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AN INVESTIGATION OF ATMOSPHERIC TURBULENCE 
BY STELLAR OBSERVATIONS* 

by 

Jack Lytle Bufton 
Goddard Space Flight Center 

CHAPTER I 
INTRODUCTION 

Light propagation through the atmosphere is a subject that until recently 
has received neither experimental nor theoretical treatment. Traditionally, the 
study of this phenomenon has been concentrated in astronomy where the collection 
of visible radiation from stellar sources is of utmost importance. Recently, the 
study of propagation has attracted the scientist interested in relating the physics 
of the atmosphere to observed effects, as well as the engineer interested in 
defining atmospheric limits for various optical systems. The advent of the laser 
and its implications for optical communication and ranging has spurred interest 
in the whole field of atmospherics. 

Although to date there is no unified treatment of atmospheric effects at 
optical frequencies, several investigations have been partially successful in pre- 
dicting experimental results. There exists a body of journal articles and research 
reports treating various propagation problems and estimating limits imposed by 


*The information presented herein was submitted as a thesis in partial fulfillment of the requirements for the 
degree of Master of Science in Physics, University of Maryland, College Park, Maryland, 1970. 
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the atmosphere. The most important areas for future investigation lie in the rela- 
tion of physical properties of the atmosphere to optical data and the determination 
of applicability and range of validity of the present theories. The purpose of this 
thesis research is to provide more information in the first area in order to shed 
some light on the second area. The approach is a study of a particular problem, 
the propagation of light from a star to an earth -based telescope receiver. Theoret- 
ical analysis is used to predict effects on starlight and experimental data is then 
used to evaluate the models and relationships behind this analysis. 

At the outset an introduction to the physics of the problem is presented. This 
is followed by sample diffraction calculations for a thin layer of atmosphere and a 
more detailed introduction into the atmospheric physics. Statistical theories of 
recent years are then reviewed and the pertinent results for starlight propagation 
assembled. Where possible, formulae of other researchers are used intact. In 
several situations the theory must be extended or modified to more closely corre- 
spond to the experimental part of the thesis. In particular an analysis of stellar 
image intensity and motion-time behavior is presented here in a new form. A de- 
vice developed at Goddard Space Flight Center to study stellar image quality in a 
telescope is used to provide quantitative data. Data is collected on stellar image 
intensity, motion, and size as observed in a 0.152-meter diameter refracting 
telescope. 

The central contribution of the thesis research is the use of data to evaluate 
parameters of a model turbulence profile. Such a profile is of general interest to 
all researchers in the vertical propagation problem regardless of their specific 
optical systems. Data on stellar image intensity and image motion -time behavior 
is combined with upper altitude wind velocity and the theoretical analysis developed 
earlier to arrive at the model parameters. A concluding section discusses applica- 
tion of the model and directions for future research. 
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CHAPTER II 


THEORETICAL BACKGROUND 
Physics of the Problem, Part I 

Optical propagation in the earth^s atmosphere is influenced by absorbtion, 
molecular and aerosol scattering, and turbulence -induced scattering. The first 
two mechanisms result in effects that can be minimized by proper choice of local- 
ity and weather conditions for observations. In addition they are relatively sta- 
tionary effects within a particular locality or in a period of several hours. The 
third mechanism, turbulence -induced scattering, results in the most serious 
optical problems from the viewpoint of this thesis. The turbulence process pro- 
duces random fluctuations in refractive index along the optical path. This in turn 
distorts the phase or shape of the optical wavefront. After propagating a distance 
through the turbulent atmosphere an optical wavefront will exhibit phase and 
ir radiance distortions. These distortions can change rapidly in magnitude over 
short periods of time and among different localities. 

First attempts at solution of optical propagation problems, ranging from 
the work of Lord Rayleigh^ to more recent investigators,^’^ suggested that the 
observed optical effects may be explained by thin layers of turbulence in the 
atmosphere. Diffraction calculations relate optical effects to the statistics of the 
layers. In reality effects from many levels in the atmosphere are important in 
producing the observed results. Without considering, for the moment, the 
turbulence process in more detail, it is instructive to examine a highly simplified 
one-layer atmospheric model. Although far from rigorously correct, the single- 
layer model does develop the trends and reveal the optical physics involved. 
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Single-Layer Model 

Consider the two-dimensional optical propagation problem, Figure 1, de- 
fined in the following manner. A unit amplitude plane wave (starlight) is incident 
upon a thin atmospheric layer that imparts a sinusoidal modulation to optical 
phase. The disturbed wave propagates to a ground-based telescope where an 
image is formed. The telescope is modeled by a one -dimensional aperture and 
thin diffraction limited lens. The quantities of interest are expressions for ir- 
radiance and phase of the wave in the aperture and image planes in terms of the 
various propagation parameters. 


WIND 

VELOCITY 


INCIDENT UNIT 
AMPLITUDE PLANE 
WAVE 


Z = 0 PHASE^ DISTURBANCE 


U(X;0) 



APERTURE PLANE 



Figure 1. Thin layer propagation problem 


The propagation parameters are as follows: 


0(x) = optical phase modulation imposed by thin layer 


<p (x)= Bcos [2n (x/d + A)] 

d = spatial wavelength of the phase pattern 
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A = A (t) = relative phase 
A= v^t 

= velocity of pattern (wind velocity) 

B = depth of phase modulation 

z = propagation distance from layer to aperture plane 
2L= aperture width 
k = optical wavelength 

k = 2tt /k . 

After passing through the phase modulation the plane wave disturbance is repre- 
sented as follows: 


U (x; 0) = e xp [ i k B co s 2 (x/d + A )] • 


( 1 ) 


This can be written in terms of an infinite series using the Jacobi-Anger expan- 
sion, 


U 



(k B) exp [ 2rr im (x/d + A + 1/4)] . 


The angular spectrum of U(x; 0) in terms of spatial frequency, f is the Fourier 
transform of U(x; 0) as given by the expression 


Ao(fx) = 



U (x; 0) exp ( - 2 77 i x) dx 


( 2 ) 


^ (k B) exp [277 im(A + 1 '4) ] S (f^ - m/d) . 
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The angular spectrum, (f*)* after propagation the distance, Z, is just (f^^) 
with phase delayed by the factor, k z Vl - ( \ , 


C^x) = A) (Q ii - (K f^)2] , 


(3) 


The desired quantity is U(x; z), the complex disturbance in the aperture plane. 
This is the inverse Fourier transform of (fx). 


1 


U(x; z) = (fj exp [ikz Vl - (\ iy]x 


(4) 


exp (2 7T i f x) d f 



(kB) exp [2rrim('A + l/4)]. 



S - m ^d) exp [ikz 


exp (2 77 i x) d 


(k B) exp [2 77 i m ( A + 1 /4)]x 

m = - 00 

X, m 


exp 


exp [2 77 i 


The last result is obtained by substitution of equation (2) in equation (4) and per- 
formance of delta function integral. With the relations 
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Ill^ 


and 


J_„ cx) = ('- ir j„ (X), 


iin X = exp (ix) - exp (- ix) 
2i 


cos X = exp (ix) + exp (- ix) 


U(x; z) can be rewritten 

U (x; z) = (k B) exp (i k z) (5) 

+ (k B) exp [ i k z 1 - ( \ /d)^] -X 

2 i sin [277-(A + 1^4 + x/'d)] 

+ (k B) exp [ikz l-4( \ ^d) X 
2 cos [4 tt(A+1/4 + x ^d) ] 


Under the assumption that the depth of phase modulation is much less than a 
wavelength, the argument of the Bessel functions is much less than one, that is, 


kB.iiLB < < 1 

X 


Thus in the series expansion of U(x; z) only terms of order (kB)^ or less should 
be retained. From the series expansions of the Bessel functions themselves, it is 
found that 
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where 


00 

Jn. ~ Jo + Jl + J2 > 

m = 0 


and 


Jo(kB)s:l-(!^J. 

Jj (k B) ^ k B ^2 , 

J 2 (k B)^i: (k B/8)2, 


to order (kB)^. An additional necessary assumption is thatX /d < < 1, i.e. , the 
wavelength of light is much, much smaller than the spatial wavelength of phase 
modulation. This is quite reasonable since the smallest dimensions of turbulent 
eddies are on the order of millimeters. Thus, 

Vl - (A. ~ 1 _ 1 (X /d')2 , 

2 

and 

Vl - 4 (\ /d)2 1 - 2 (\ . 


Employing the above approximations, equation (5) results in 


U (x ; z) ~ exp (i k z) 1 - 


k2 B2 


( 6 ) 


+ k B i cos [2 77- ( A + X M)] exp 


i k z 


(X /d)' 


k2 B2 


cos [4 77 (A + X ^d)] exp [- i 2 k z (X/d)^] 
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Ir radiance in the aperture plane is now formed, as follows: 


I (x; z) = U* (x; z) U (x; z) 


k2 

I (x; z) = 1 - — - — + 2 k B cos [ 


2 7T (A + x/d^llsin ^ 


( 7 ) 


k2 B2 


COS [47T(A + x/d)]cos [2kz(\/d)2] 


k2fi2 k2B2 .. .. 

4- + COS L477 (A + X /d) ] . 

2 2 


To proceed further a small-angle requirement is imposed. Thus, 


k 2 ( \ ^d) 2 < 77 / 1 0 , 
sin [ k z ( ^d) 2] ~ k z ( \ /d) ^ , 

and 

cos [k z ( \ ^d)^] ^ 1 . 


For z = 10 meters and \ = 0.5 micron this forces d to be greater than one milli- 
meter. For z = 10 kilometers and ^ = 0.5 micron d must be greater than 10 centi- 
meters. Thus the spatial wavelength, d, which corresponds to turbulence dimen- 
sions, must be larger than a certain minimum size which is near a millimeter at 
the earth’s surface and rises to approximately a meter at the edge of the atmos- 
phere. The data or estimates available on minimum turbulence sizes as a 

function of altitude suggest that these requirements are reasonable. Thus I (x;z) 
reduces to 
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(8) 


I (x; z) ^ 1 + cos [2 t7(A + x/cO]. 

d2 

The original sinusoidal modulation of phase has resulted in a sinusoidal irradi- 
ance modulation in the aperture. The depth of this modulation is 


4 rr^ B z 
d2 


and thus increases with growth of z and B and decreases for larger pattern sizes 
at a given z. 


Total intensity in the image is found by integration over the aperture as 
follows: 


TOT 


I (x; z) d X 


T ^ ^ A \ 

2 L -f cos (2 t7A) SJLn 


(9) 


2 7T L 


Restriction to the dimensions of the aperture introduces the additional term 


. 2 77 L 

sin 

d 

that acts as an aperture -averaging term. Changes in 2L, the aperture size, 
are not as important in effecting where 2/d >> 1 as when L/d ~ 1. > 

like I(x; z), is actually a function of time through the dependence of A on time, 

A = v^ t = A (t^ . 
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Thus, the frequency of image intensity variations (scintillation) is directly pro- 
portional to wind velocity at the particular altitude, z, of the phase modulation 
layer. 




This same type of analysis can be applied to arrive at an expression for opti- 
cal phase in the aperture plane. Together with irradiance information, a trans- 
formation to the image plane can then be made with the Fraunhofer dif- 
fraction approximation. The interesting quantities of image size, image motion, 
and their respective time dependences can be computed in terms of the initial 
propagation parameters. The interested reader is referred to the work of 
Ratcliffe^ and Ramsay"* for detailed treatments in this area. This analysis will 
not be presented here. Instead a presumably more accurate analysis based 
on current theories of propagation and the actual structure of turbulence 
is presented in Chapter III. Although inaccurate for specific results, 
the foregoing analysis is not without merit. As irradiance behavior in the 
aperture plane is again discussed in Chapter El in light of current theories, 
many of the concepts developed from the simple one-layer model will be 
recognized. 

Physics of the Problem, Part II 

Atmospheric optical index of refraction, N, is not sensibly affected by 
pressure independent of its change in density. In the visible region there is a weak 
dependence on wavelength. The following equation expresses these dependencies 
and includes empirical constants established^ for the visible range of the 
spectrum: 


N = 1.0 + (^288 + L— 10"® , 

^0 V y 


( 10 ) 
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where 


p = atmospheric density, 

= (standard temperature and pressure) density of dry air, 

and 

\ = optical wavelength in microns. 

Atmospheric density at a given point undergoes rapid, small fluctuations around 
a relatively steady mean value. This is a result of turbulence mixing the air. The 
causal phenomenon for turbulence may be thermal or mechanical. Detailed analy- 
sis of these mechanisms is left to the meteorologist or atmospheric physicist. 
From the optical standpoint the important points are that turbulence results in 
quite serious optical problems and that its structure can be described with some 
success. 

In typical turbulence, generation energy is input to the atmosphere in large- 
scale structures, with dimensions of meters or larger. Energy is then trans- 
ferred to smaller and smaller turbulent eddies with little energy loss until finally 
it is dissipated into heat in the smallest eddies (dimensions near one millimeter 
at the earth^s surface). Turbulence characteristics between these ill-defined 
outer, Lq, and inner, scales which make up the "inertial subrange" can be 
treated in a reasonable fashion. 

It is evident from the random nature of turbulence -induced fluctuations that a 
statistical approach is called for. The primary statistic of the random field is 
the structure function, defined as the mean square difference of a quality of the 
field for two points separated a distance r . If the field can be considered homo- 
geneous and isotropic the structure function depends only on the magnitude, I r | , 
of T . This will be denoted simply by r. 

Analysis due to Obukhov^ and Kolmogorov ^ concludes, primarily on the basis 
of dimensionality arguments , that the structure function for velocity fluctuations 
in the inertial subrange can be written 
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(r) = <^[v (7j) - V 

= c (e r)2/3 , -^0 \ ’ 

where 

D^ (r) = velocity structure function , 

r = separation of measurement points 1 and 2, 

and 

e = energy dissipation rate per unit volume. 

Associated with (r) is its Fourier transform, F^ (k), the velocity spectrum. 

In three dimensions F^ ( k ) is proportional , where ? is a spatial frequency 

vector. This formulation has been verified experimentally®-^ by measurements 
of small-scale velocity and temperature fluctuations in atmospheric surface layers 
and shown to be generally independent of the mechanism used to generate the 
turbulence. However, density, not velocity, is the important parameter from the 
optical viewpoint. One step away from density fluctuations and, thus, refractive 
is the temperature structure function, D^(r). Obukhov concluded that (r) can 
be written in a form similar to (r). Thus, 

(O = C2 r 2 '3 , -Cq << r << (11) 

where 

= temperature structure constant. 

Temperature and density statistics are directly related since the atmosphere will 
not support local pressure gradients. Warm segments of air are less dense than 
similar segments of cold air at the same altitude. Typical temperature fluctua- 
tions are less than a few tenths of a degree centigrade. The resulting density 
changes are of the order of hundreds of a mg/cm^ and refractive index fluctuates 
by one part in 10“^. The additive effect of many small changes of N (p;\) along 

the optical path can result in serious wavefront distortion at the receiver. 
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Mathematically the connection is as follows. The variation, A p, in 
density (at constant pressure) resulting from a small variation in temperature, 
A T, is 


A p = 


P 

RAT 


from the ideal gas law. Further, since 


lim A p _ p _ p 

AT-^O AT R T^ ^ 

then, 

A p ^ A T 

H ^ 


From equation (10) the effect on refractive index of a small change in density, 
Ap, can be observed. Thus, 


AN (p; A) = 1 . 0 + Af f 288 + ] x lO"® . 

Po V A2 / 

Statistics of A N and AT are thus connected by the equation 


AN (p; \) = 1. 0 - — (^288 + iAlA lO'^. (12) 

Po T V y 

Structure functions relating N and T go as the square of A N and AT, Thus the 
two structure constants are related in the following manner: 





288 + 


1 . 34 


X 


2 


10"6 



(13) 
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where 




D2 (r) . C2 r2/3^ ^^,«r«Lo. (14) 

These equations are the important results from turbulence theory. They relate 
the statistics of refractive index to a strength of turbulence factor, , and in turn 
to the strength of temperature fluctuations, - 

In the foregoing analysis isotropy has been assumed. For separations, r, 
well within the inertial subrange and at the same altitude this is not a bad as- 
sumption. Strength of turbulence, , takes on different values as a function of alti- 
tude due to changes in absolute density, absolute temperature, and the strength of 
temperature fluctuations. For vertical propagation problems, altitude dependence 
must be taken into account. This is accomplished by assuming local isotropy at a 
given altitude, retaining the ^' 2 / 3 " law for (r), and considering to be (h), 
a smoothly varying function of altitude, h, as follows: 

(r: h) = C2(h) r2/3^ « r « Lg . (15) 

Similarly sizes of inner and outer scales vary with altitude. The inner scale is 
thought to attain its smallest value near the ground and increase approximately 
as . The outer scale is less well understood but may be approximated by 

Lq 

where 

h = height in meters . 

These rough values for and were obtained from Hufnagel^^ and Fried 
Various estimates and data are also available on the C,^ (h) profile. These are 
presented as an introduction to Chapter V. 
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Finally, it must be noted thatthere exist turbulence regimes (particularly at 
night near the earth’s surface under conditions of low wind) where the existence 
of an inertial subrange and the applicability of the "2/3” law are questionable. In 
these regimes there is insiiEficient energy input to have well-developed turbulence. 
For altitudes above the surface layer the existence of high Reynolds number flow 
and the associated well -developed turbulence argue for the existence of inertial 
subrange. 

Introduction to Current Theories 

Expression of optical effects in terms of the statistics of refractive index 
has been the subject of many current theoretical treatments. Since the earlier 
work in diffraction from single layers, researchers have employed geometrical 
optics and approximate solutions to the wave equation. In geometrical optics 
the usual starting point is integration of the "ray equation" as in Chernov^^, 
Beckmann^^ , and Hodara^"^. Along the way, assumptions of no diffraction and 
small -angle scattering are made. Range of validity of this approach is restricted 
to the turbulence pattern near field. Since inner scale of turbulence, is the 
smallest dimension of the pattern, the range constraint can be expressed 

0 v/'\ z , 

where 

z = range of validity. 

This is quite a severe constraint, of the order of a few meters near the ground. 
Recently a few authors, notably Taylor^ have indicated that geometrical optics 
may apply for longer ranges, at least for predicting optical phase statistics. If 
a Gaussian form for refractive index structure function is assumed, the theory 
lends itself to use of geometrical optics. This is where considerable effort has 
been placed in the past. Here the important parameter becomes a correlation 
length, , intermediate in length between 'i ^ and L^. Much evidence cited 
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previously suggests that the Obukhov-Kolmogorov (OK) theory is more applicable 
especially near the earth’s surface. In addition, geometric optics results for the 
Gaussian case have been shown to be somewhat in error by Kerr^®. Application 
of geometrical optics and the OK turbulence spectrum to stellar observations is 
presented by Rieger and discussed further by Young 

The most successful treatments of optical propagation through turbulence have 
been achieved by approximate solutions of the scalar wave equation 


V2 u (?) + k2 (?) u (?) = 0. (16) 

For a plane wave incident on the turbulence, as in the case of starlight at the top 
of the earth’s atmosphere, the complex field distribution is written 

u (?) - A (?) exp [i k 0 (?)] 


or 


u (?) 3 exp [in A (?) + i k 0 (?)] . 


(17) 


As this disturbance propagates through the atmosphere it is randomly perturbed 
in a multiplicative manner. The end result is a product, not a sum, of the many 
intermediate actions of randomly fluctuating refractive index. Terms in the 
exponent of equation (17) are effected in an additive manner. As a consequence 
of the central limit theorem these terms should have a Gaussian distribution. 

Thus the conclusion follows that log-amplitude and phase obey Gaussian statistics. 


Of the wave equation treatments use of the so-called ’’Rytov approximation” 
has achieved preeminence. The basis of this approximation is the neglect 
of a small term in the wave equation after the expression (17) is substituted in. 
This transforms it into a linear equation which can be solved for the desired 
optical statistics. Principal investigators of this method are Tatarski^^ and 
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Fried There is extensive literature devoted to these mathematical efforts. 

For a full mathematical treatment the interested reader is referred to Tatarski’s 
"source book". While arguments are still raging over the validity of this method, 
it has been verified experimentally 21^22 least for scintillation over ranges be- 
yond the geometrical optics limit. Both American and Russian experimenters 
have, however, reported saturation effects for long paths and/or high turbulence 
conditions near the ground. For vertical propagation no such saturation effects 
have been observed and it is felt that the theory holds well. As with other optical 
propagation theories direct confirmation eludes the observer since the theory is 
really a two-step process. First the optical statistics of interest must be related 
to refractive index statistics of the field of turbulence, such as in the Rytov ap- 
proach. Then the turbulence process must be adequately described. An error 
in either or both would have an adverse effect on the prediction of experimental 
data. This is especially true for the vertical case where the strength of turbu- 
lence changes with altitude. 

Another school of investigators, notably Beran^^ and HufnageP'^, while work- 
ing from the wave equation, have modified it to a differential equation in statistics. 
Averages are taken before the equation is solved. The end result is a mutual 
coherence function for the atmosphere. This can be directly expressed in terms 
of Tatarsksi^s wave structure function. Thus, at least in the area of average wave- 
front degradation, the two principal methods achieve the same end result. Once 
again prediction of experimental data depends on the form of the wave structure 
function. This in turn depends on the representation of turbulence. 

In the following chapter results of theories applicable to stellar observations 
will be discussed. In several cases the present theory must be extended or 
modified. 
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CHAPTER III 


ANALYSIS OF STELLAR IMAGE EFFECTS 


Image Intensity 

Total image intensity results from the irradiance pattern of starlight sampled 
by the telescope aperture. Irradiance varies with both position and time in the 
aperture plane. As indicated in Chapter II the quantity log-amplitude should obey 
Gaussian or normal statistics. Thus irradiance should have a log-normal be- 
havior. Log-amplitude is defined in te^ms of irradiance, by the equation 


I (x, 



t) 


( 18 ) 


where 

^ (x, t) = instantaneous log-amplitude at position x in the aperture plane, 
I (x, t) = corresponding irradiance at x , 


and 


Iq = ensemble average of I (x, t) 


= <i (X, t)> . 


Variance of log-amplitude, cr 
average log-amplitude 


2 


has been shown^^ to equal the negative of 


O'! = - (X, t)> 

under the assumptions of Gaussian statistics for log-amplitude and conservation of 
energy. A theoretical expression for cr^ has been proposed by Tatarski^® on the 
basis of the Rytov approximation and Kolomogov turbulence spectrum. Thus, 
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(19) 


^CO 


cr^ = .56 sec 011/6 


(h) h®/® dh 


where 


k 


= optical wave number = 


2 7T 

~ 


and 


6 = stellar zenith angle , 

(h) = refractive index structure constant as a function of altitude, h. 


This formulation is also conditioned on an infinite plane-wave source (starlight) 
and a point detector. A telescope with small aperture such that instantaneous 
irradiance is constant over the aperture can be considered a point detector. The 
predicted wavelength dependence contained in k^^^ is such that is reduced 
by 20 percent in passing from A. = 0,45)Lito\ = 0. 55/^. While cr^ is reduced in magnitude 
at the higher wavelength, analysis has shown that spectral correlation of 
remains high for this separation of wavelengths. 


Secant of zenith angle is simply the air mass through which starlight must 
pass before reaching the telescope. The above 11/6 power dependence has been 
shown to break down for large zenith angles when detection is not made under 
monochromatic conditions. A recent study has attributed this to the increasing 
role of chromatic dispersion for large 6 . Under such conditions different wave- 
lengths of starlight suffer varying degrees of refraction in the atmosphere so that 
they follow different paths through spatially uncorrelated turbulence. The end 
result is a saturation -like effect for at large 0, Dependence on strength 
of turbulence, (h), involves the 5/6 moment of altitude. Thus altitudes where 



(h) dh 
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has its largest values are most important in determining The significance 
of this important relationship will be discussed in Chapter V, 


A photomultiplier tube or other square -law device is normally used in field 
work to obtain data on the stellar image. Electric current out of this device is 
proportional to irradiance, not amplitude. Thus it is sometimes more convenient 
to handle the data in terms of irradiance statistics. The connection with log- 
amplitude statistics is easily made. Let the ratio of rms image intensity fluc- 
tions to their average value be the fundamental measure of irradiance statistics. 
This ratio, here called CIV, is the coefficient of irradiance variation. Thus, 


where 


CIV 



0-2 = irradiance variance 


= <(I - Io)^> • 


Noting that 


and, thus, 


0-2 can be written 



I = Ioexp(2-e), 


( 20 ) 


( 21 ) 


°~i = exp 

= Iq ^[exp (2-&) - -(exp (2^))> ]2)> , 

where 

<exp (2^)> = exp [2<-e> +2 )2> ] 

= exp [2 + 2 o-^] 

= 1 , since ^ 
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The last result follows from the Gaussian behavior of •i. Likewise, 


and 


<(exp (4't)'^ = exp [4 + 8 cr|] 

= exp (4 o-^) 

- ^0 <(l^exp (2'C) - 
= <(exp (4't) - 2 exp (2^) + 1^ 
= [exp (4o-|) - 1 ] . 


The ratio 0-2/12 is just CIV^, Thus, 


CIV2 = exp (4o-^) - 1 . 


( 22 ) 


In actual experimental situations where equipment dynamic range or small de- 
partures from Gaussian behavior must be accounted for, the use of this equation 
is questionable 22 . For values of o-^ <0.05 where there is little difference between 
the normal and log-normal distributions the equation may be used safely, however. 
Furthermore, for these smaller values an approximation to equation (22) can be 
employed as follows: 

CIV2 = exp (4cr^) - 1 
CIV2 + 1 = exp (4a|) 

In (CIV2 + 1) = 4cj| 


For 


Therefore, 


= CIV^ 



CIV^ < 1, 

In (CIV2 + 1) :^CIV2 . 
Civ2 4cr| . 


(23) 


Experimental studies of starlight almost invariably employ telescopes of 
aperture dimension larger than several tens of centimeters. These optical sys- 
tems cannot be considered point detectors. In the experimental phase of this 
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thesis the aperture is circular and unobstructed with a diameter of 15.2 centi- 
meters. This is intermediate between the point detector and a typical astronomi- 
cal telescope of a meter or more in diameter. As diameter increases the 
strength of image intensity fluctuations, measured either by or CIV, decreases, 
slowly at first and then more rapidly as the spatial correlation distance of ir- 
radiance in the aperture plane is exceeded. In addition, the frequency content 
and zenith -angle dependence of the fluctuations are altered. 


The interesting statistic now becomes covariance rather than variance. 
Covariance of log-amplitude, cr|(p) , and covariance of irradiance, 0-2 (^) , are 
functions of aperture plane dimension, p , where 


and 


0 £ p $ D 

D = aperture diameter. 


Explicit covariance results for the starlight propagation problem, based on the 
Rytov approximation and Obukhov-Kolmogorov turbulence spectrum, are limited 
to those of Fried He has achieved results only by assuming a specific form for 
C^ (h) as a function of altitude. More general results corresponding to equation 
(19) are unavailable. The work of Fried and others, notably Tatarsksi^^ and 
Yura and Lutomirski^® still provides a solid basis for understanding the problem. 
They all suggest that the aperture averaging effect can be described by a 
functional dependence involving diameter, D, correlation distance, , and 
strength of fluctuations for the point detector case. Thus the variances for log- 
amplitude and irradiance measured with a telescope of diameter, D , can be 
written 




= A 


■I 


./\z 




( 0 ) 


( 0 ) 
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or 



— A 

D . (0) 

o to 

1 ~ 

D 

/ATz T2 

- -*-o J 


( 0 ) 

12 

■^0 


[CIV2]j^ = Aj 


L v'XT 


CIV- 


CIV2 . 


( 24 ) 


The correlation distanceX^ z , is expressed as sec 6 for stellar observa- 

tions at zenith angle, 6 . Altitude, can be considered a scale height for the 
atmosphere, the altitude below which most of the distortion-causing turbulence 
lies. The factor z arises from the relation of two-dimensional spectral density 
of amplitude, , to three-dimensional spectral density of the structure func- 

tion of refractive index, (x) . As given by Tatarski^^ this relation is 


(K)^^277k2 


^co 

Cn (z) 



dz 


(25) 


for vertical propagation, subject to the condition that 


^^«Az«Lo. 

F^(/<) has its maximum near k = l/A z , The two-dimensional Fourier Bessel 

transform of F. {k) is log-amplitude covariance, (p) , in the telescope aperture 

A a, 

plane. Thus, 


<v| (p) = 2 tt 


Jq(/< p) (k) k dk 


(26) 


p 00 roo 

= 4772k2 C^(z) 

A •'a 


(k p)y 


2 

"2ir 


<dK d z , 


where 

< I ) f ^ (k) = 0.033k~^^^\ 
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A peak near k - I//X 2 in the spatial frequency domain means will have 

a correlation distance on the order v/AT. From Fresnel diffraction theory ^/kz 
will be recognized as the radius of first zone. Thus the general statement can 
be made that turbulence of dimensions v/\z located a distance z from the receiv- 
ing plane will have the largest effect on o-^ (p) . 

The quantity can be inferred from a direct measurement of (p), 
although such data is difficult to obtain. A second method based on aperture 
dependence of CIV^ will be used instead. Data of this form is available as 
reported in Chapter IV. Use is made of equation (24) for CIV^ at aperture 

diameter D to obtain 


[CIV2]^ A, 


^/}^z 


CIV^ 


CIV2 . 


(27) 


The aperture averaging factor Aj is roughly determined by Tatarskiis equation 
(13.28) for a quantity he labels G. Numerical integration of G as CIV^ approaches 
zero is presented in his Figure 36. For values of D//Vz between 0.2 and 3.0, 

G can accurately be modeled by an exponential. Thus, 

G = 1.20 e" d//7Tz 

(28) 


and 


= 1.20 Civ2 . 

Experimental data for CIV^ as a function of aperture, D, is then also fit with 
an exponential by the expression 


[CIV2]j^ . <29) 

Equating these two expressions for [CIV^lp yields a solution for in terms 
of the data -fit parameter, B. B has dimensions of length"^ , such that 


/n = - 1.33/B. 
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This becomes 


v'X.H 


0 ” 


- 1.33 

B v/sec 6 


(30) 


when allowance for stellar zenith angle, i9, is made. 

Since Aj also depends on 6 as well as D through the term Hq sec (9 , 
the zenith angle dependence of [CIV^] ^ will change for various ranges of the 
ratio T>/ V K Hq sec 6, Tatarski points out that for small values of 0, [CIV^Jq 
can be approximated by a sec^ & , where /3 is a power dependent on aperture 
size. As D increases so should from the point detector value of 11/6 to a 
value of 3.0 when D » sec <9 , 


The effect of exceeding the correlation distance has been measured by 
Protheroe^^. He obtained a 2.5 times reduction in CIV in going from a 2,54-cm 
diameter aperture to al5.2-cm diameter. This corresponds to a reduction in 
(p) of 6.25. For the 15. 2 -cm -diameter aperture he observed the following 
dependence on zenith angle : 


and 


CIV « seci-2 e 


« sec^'^ 6 . 


(31) 


The calculated correlation distance was about 8 to 10 centimeters, Protheroe's 
actual CIV values must be corrected for this experimental technique in order to 
arrive at CIV for zenith viewing. When this is done the following average values 
can be presented; these in turn predict average values. All measurements 




were taken with a 15. 2 -cm aperture. 



CIV 

CT^ 



;L 

Summer 

0.14 

0.0049 

Winter 

0.16 

0.0064 
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The frequency spectrum of irradiance is determined by the same factors of 
correlation distance, wavelength, zenith angle, and aperture size. In addition, 
the velocity spectrum of upper altitude winds becomes very important. Increases 
in zenith angle or aperture size are reflected in a decrease in magnitude of high 
frequency components. These effects are observed in spectral densities of stellar 
image intensity and log- amplitude. 


The physical basis of spectral -density aperture dependence can be understood 
by considering the structure of irradiance in the telescope aperture plane. The 
structure consists of a continuous size spectrum of bright and dark patches 
changing with respect to time and sweeping across the aperture. As the aperture 
becomes much larger than the correlation distance the small patches, which have 
dimensions of the order of the correlation distance and are responsible for the 
high frequency components, produce less of an effect due to aperture averaging. 
Larger particles and, therefore, lower frequency components do not suffer as 
much averaging. The rate at which the pattern or parts of the pattern are swept 
across the aperture is determined by upper altitude wind profile. In the middle 
latitudes wind velocity is strongly peaked in the 8- to 14- kilometer 
region of the tropopause. This is not surprising since turbulence generation is 
readily associated with high flow velocities. Typical velocities at the tropopause 
are from 20 to 100 meters per second and the flow is often referred to as a jet 
stream. On the basis of the above physical considerations, a proportionality to 

irradiance frequency can be asserted as follows: 


f cc 


(32) 


/\Hq sec d 


where 


V = component of wind velocity normal to the propagation path. 
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Theoretical analysis of irradiance frequency behavior best proceeds in the 
manner revealed by Tatarski. His equation (13.38) expresses spectral density 
of image intensity in terms of an integral over the two-dimensional spatial 
spectrum of irradiance and a filter function for a circular unobstructed aperture 



2 Jj (d/ 2 V /<2 + 477^ fVv2) 

(D'2 VK2 +4 7t 2 f2/V2) 


2 

dK, (33) 


where 


and 


= component of wind velocity normal to the propagation path , 
K = spatial frequency in aperture plane , 
f = temporal frequency, 

D = aperture diameter , 

Iq = average irradiance, 

Fj (k) = two-dimensional spatial spectrum of irradiance. 


To arrive at this point Tatarski has employed the "frozen in" theory of turbulence 
structure. The physical basis of the theory is the assumption that wind transports 
turbulent structure across the propagation path before the structure can itself 
change. In meteorology this would be equivalent to the assumption of Taylor^s 
hypothesis.^ The result is that irradiance, I(x , t + r ), at point, x , and time, t , 
plus a delay, r, can be simply connected with I(x, t ) displaced by a distance 
Vj^T by the expression 

I (x, t + r) = I (x - r, t). (34) 

Autocorrelation of irradiance, Rj (r), can then be expressed with this result as 
follows: 


Rj (r) 


<([I (x, t + r) - Iq] [I (x, t) - Iq])> 

<I (X, t)2> 


(35) 
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Equation (13.38) then results when' correction for circular aperture is made and 
K is shifted by an amount 2 tt f/v„. 


Taylor’s hypothesis has been verified experimentally by several researchers. 
For the particular case of vertical path the reader is referred to the work of 
Gossard^^ . A natural question at this point would concern the range of validity 
of Taylor’s h3q5othesis. Tatars ki has shown that the condition for validity can 
be expressed as 


where 

and 


Lq = outer scale of turbulence, 

z = distance from detector to the turbulence generating region. 


That is, the dimensions of turbulence structure that effect light propagation should 
be much less than the dimension of the flow as a whole. This is simply one side 
of the condition for operating within the inertial subrange. 


Tararski development of (f ) beyond equation (13.38) assumes a constant 
value for strength of turbulence and wind velocity, as might reasonably be ex- 
pected for a horizontal path. Observations of starlight require a knowledge of 

the interaction of all regions in the vertical path. These regions contain different 
2 

values of and and must be weighted according to their distance from the 
receiver. This can be accomplished mathematically by obtaining the vertical 
propagation relation for Fj(/<), inserting the expression into equation (33), and 
interchanging the order of integration so that the k integration can be performed 
first. Fj (k) is obtained from F^(/<), the two-dimensional spectral density of 
amplitude in equation (25) , by noting the relation 


F, (K) (k) . 
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Thus, 


Fj (/<) ^8 t 7I2 k2 C2 (z) (/<) sin 

•^0 

J ^oo 

I 

r, 


2 Z 


2k 


dz 


.264 7T J2 k2 


(z) sin^ ( — -] dz. 

2k 


( 37 ) 


This is the vertical propagation relation for Fj (k). 

Inserting this approximation into equation (33) yields 


W (f) Ji; 2.11 T? k^ 


^ 00 ^00 

1 i 


(z) 


2 4t7^ 

+ 


- 11/6 


(/<2 + 4t7-2 ^ 


= 2. 11 77^ k 


2k 


2 w2 


dz 


(z) 


2Jj (D/2 V /<2 + 4 tt 2 fVV^ (z ; 


X 


12 


r eg ( z ) r 
Jo \ (^) Jo 


D/2 y + 4 77^ 


[/<2 4- 4 77^ X 


s in 


-1- 4 77^ ^ 


2Ji (D/2 //<■" + 4 t 7^ fVv2 (z)) 

2k 

L J 


D/2 + 4 77^ f ^/V ^ (z ) 


d K dz 


poo Q 2 / \ 

2.11 77^ k^ — F {f, (z), D, k} d z. 

•'o ''n 


(z) 


(38) 


This result shows that irradiance frequency spectrum is expressable in terms of 
an integral over the turbulence and velocity profiles and a non-linear function F 
{f , V^(z), D, k}. The specific form of F depends on the interaction of wind 
velocity with aperture shape and the spatial spectrum of turbulence fluctuations . 
In Chapter V this result is used to construct a model forC^ (h). 
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Image Motion 

Fluctuations in the position of image center of gravity result from atmos- 
pherically induced changes in the phase or shape of the incoming wavefront. 
Most changes can be attributed to random linear phase shifts or tilts in the 
wavefront.^® The other principal wavefront degradation is quadratic in form 
and leads to focusing errors, thus contributing to blurring of the image. For 
a point detector, wavefront tilt will result in displacement of the entire image. 
Tatarski^^ has related wavefront tilts to phase structure fimction by examining 
a two-point correlation (interferometer) experiment, Figure 2. 



ANGULAR DISPLACEMENT IN 
IMAGE PLANE 

Figure 2. Interferometer case 

Consider two point apertures, A and B, separated by a distance, d. Incoming 
starlight represented by a wavefront tilted with respect to d is brought to a focus 
after passing through the point apertures. Angle of tilt, A a radians, is equal to 
the angular displacement of the diffraction pattern in the image plane. Optical 
phase shift, Ac^ , between A and B, can be expressed in terms of this tilt angle by 
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A0 = (jbg-0^=kdAa, 


(39) 


where 



Since the phase structure function, (d), is defined as the mean square 
difference of phase between two points , 

<(^B - (d) • 

Then, the above relation for tilt angle yields the expression 

<(Aa)2> = . (41) 

^ ^ k2 d2 k2 d2 

The angular brackets once again denote an ensemble average. Tatarski^® em- 
ploys his expression for phase structure function based on the Rytov approximation 
and the (OK) turbulence spectrum as follows: 


(d) = Kk2 d®/3 



(z) d z • 


(42) 


The value of K depends on the relationship of d to the correlation of irradiance 
fluctuations , so that 

f2.91 d > /~ 

Li. 46 d < /X.Z 


For stellar observations distance, z, along the propagation path becomes 


where 

and 


0 = zenith angle 
h - altitude. 


z = h s ec , 


Then, 


(d) = Kk^ d^^^ sec 9 


C2 (h) dh 
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0 


(43) 





Substitution in equation (41) leads to the following expression for mean square 
angular fluctuations: 

^00 

(Aa)2 = sec & (h) dh . (44) 

This is just the variance of the center of gravity of the diffraction pattern about 
the optic axis. If the stellar observation is made with a telescope instead of an inter- 
ferometer, wavefront tilt from all points in the unobstructed aperture will contribute 
to image motion. For large apertures, wavefront tilt in different regions of the 
aperture will be uncorrelated and lead to a reduction in image motion. In addi- 
tion, a non-uniform irradiance distribution will give unequal weights to tilts in 
the different regions of the aperture. As the aperture dimensions exceed the 
irradiance correlation distance this effect will diminish in importance since the 
average irradiance pattern predominates over the fluctuating part. The resultant 
image for a reasonably large aperture telescope (diameter greater than a few 
centimeters) will be blurred and exhibit a center of gravity displacement. 

In generalizing to such a case from the interferometer result above, Tatarski 
asserts that the maximum aperture dimension plays the role of the point aperture 
separation, d, with the only other change being in the value of proportionality 
constant, K. His intuitive feeling is apparently correct, as will be shown in the 
following analysis due to Hufnagel^^. The analysis will be extended here 
to include a numerical prediction of image motion variance. 

The instantaneous position of image center of gravity, r (t ), is just the first 
moment of the image intensity pattern 
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X I (x, t) dx 


f (t) = 



I (x, t) dx 


( 45 ) 


where x is a two-dimensional vector in the image plane centered on the telescope 
optical axis. Starlight incident on the telescope is represented by the complex 
field distribution, u (v , t ) , in terms of v, a two-dimensional vector in the aperture 
plane. Neglected here is the optical frequency time behavior of the incoming 
wavefront. The imaging process from aperture to image plane as performed 
by the telescope optics is described adequately by the Fraunhofer diffraction 
formula 


U(x, t) 1= C I u (v, t) exp v * x ) dv 

A. F 


( 46 ) 


where 


I(x, t) = U*(x, t) U(x, t). 


The normalization constant, C, is determined by the telescope geometry and incident 
light energy. For a lens of clear circular aperture and normalization to unit 
energy, 


and 


I(x, t) 


1 


where 

D = aperture diameter, 
F == lens focal length, 

and 

A. = wavelength of light. 


C =- 


\DF V^TT 
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As demonstrated in Appendix A and discussed in the next section on image 
size, the optical transfer function in terms of image plane coordinates can be 
written 


T(f, t) = I (x, t) exp (- 27 t i f • x) dx, (47) 

v'-oo 

where 

f = two-dimensional spatial frequency vector . 

Image center-of -gravity can be expressed in terms of T(f , t) by taking the first 
partial with respect to ? and evaluation at f =0, as follows: 


Since 


^T(f, t)_ 


= - 2tt : 


Bf 


7(t) 



t) exp (~ 27Tif*x) dx. 


t) dx , 


then 


7(t) 


i 3T(f, t) 

277 


3f 


f = 0 


(48) 


In order to express r(t) in terms of irradiance and phase in the aperture plane, 
it is first necessary to express T(7, t) in these quantities. The Fraunhofer 
formula, equation (46), can be inserted in equation (47) and the mathematics 
performed to achieve the desired result, as follows: 


I(x, t) = U*(x, t)U(x, t) 


= C 


n oo 

u*(w, t) u(v, 

OP 


t ) exp 


277i 

Tf 


i (w - v) • X 


dw d V . 
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Hence, 


T (f, t) = u(w, t) u(v, t)exp 


2 7T i 


KF 


(w-V“\Ff) *x 


dw dv dx 


u* (w, t) u (v, t) exp [ • ] dx dw dv 




or 


= C2 


-00 _00 


U (W, t) U (V, t) S 


w-(v+X.Ff) 


dw dv, 


T (f, t) 


7/D^ 


00 

I “ 


(v + A.F f ) u (v) dv , 


(49) 


The complex disturbance in the aperture plane u (v) is now expressed in terms of 
amplitude and phase by 


u (v, t) = v/'I (v, t) exp [i (v, t)], 

where 

0 (v, t) = optical phase in radians 

and 


u* (v + X.FF, t) = v'l (v + \Ff, t) exp [- i c^(v + A.Ff, t)]. 
The partial of T(f , t ) with respect to f can be performed as follows: 


BT (f, t) 4 

B f 771)2 


u (v, t) — =r u* (v + X.F f , t) dv 
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^IL “ r ‘ 


- i \F /l (v + \Ff, t) exp [- i (v + X.F f, t)] (v + t) 

3 (v + \ F f ) 


+ exp [-i0(v+X.Ff, t)] —=r v^I(v+\Ff, t)[ dv, 

3 f 


5T (f, t) 
3 f 


7=0 


4 

77 




(50) 


+ — ^ v" I (v, t ) I (v + X. F f , t ) 7 d V . 

3 f J 


The final form relating image center of gravity to irradiance and phase is 
obtained by noting that the integral over all v of the second term vanishes since 
I(v, t) has a zero value outside the aperture and by substituting the remaining 
term in equation (48) 


r (t) = I (v, t) <p (v, t) dv . (^1) 

ttD^ J_oo 

Image position thus depends on the derivative of optical phase, which is the mean- 
ing of wavefront tilt, weighted by the irradiance pattern. At this point the variable, 
t, will be suppressed for ease of notation and the variance of F calculated. 

The distribution of tilts is symmetrical about zero tilt, thus the average 
value of image motion, 7 , will be zero. The variance of image motion ex- 

pressed in radians squared, through division by F^ , is 
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00 


(52) 


<r-7>=< 


^ 6 k 2 p r 

77 ^ D J=00 


T.-vt/~^X ^0(v) 3 0(w) 

I(v)I(w) ' dv dw 




3 w 


where 

r is now measured in radians . 


The angular brackets once again denote an ensemble average. In practice the 
ensemble is formed by the time behavior of r (t). Convincing evidence exists, 
both in the literature^^ and in the author's experimental work, that irradiance 
and image motion are sample functions of independent stochastic processes. 

In that event 


n oo / 

00 \ 3 V 


(v) 3 0 


^3 V 3 w . 


3 Wy 


If in addition, the processes are considered stationary, w = v + y and (Jw = dy. 
Thus, 

16 k 2 


<V-7> - 


r f I(v)I(v+y) B0(v+yK H-Hv 

J_oc*^-co \3v 3y/ 


Because of stationarity, the expression 


3 0 (v) 3 0 (v + yr 


is invariant under translation to 



The expression for is of the form 



'To(y) Rv<?b(y) "^y > 


(53) 


where 

Ryjy) = autocorrelation of tilt, V0 

and 

TQ(y) = autocorrelation of irradiance. 

The additional assumption of isotropy allows the autocorrelations to be expressed 
in terms of 1 y |. Tilt autocorrelation is simply related to phase autocorrelation 
through the vector identity 

Rv,^(lyl) = R.^(lyl)- (54) 

Tq( |y| ) can be modeled fairly accurately by the assumption of constant irradiance 
over the aperture. The integral 

I (v) I(v + y) dv 

reduces to a result for the overlap of two circles of diameter, D. This is equiva- 
lent to the diffraction-limited modulation transfer function for the optical system. 



Phase autocorrelation is obtained from 

D^(|y|) = 2[R^(0) -R^(|y|)] 

and the Tatar ski expression for the phase structure function, 


To(lyl) 
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D,(|y|) = 2.91 k2 |y|®^^sec0 I C2(h)dh, 


J 


SO that 


^^R^(ly|) = - yI 2.91 k2 sec 0 I C2(h)dh. (56) 


^ CO 


Insertion of equations (55) and (56) into equation (53), conversion to polar coordi- 
nates, and simplification leads to the following expression: 


= j |cos“i (5)“ 

sec^f C^(h)dh. 

^ n 


2/3 


y dyx 


(57) 


The first integral can be transformed into 

r [cos”^ p - p(l - dp. 

•^0 

This has been evaluated numerically on a digital computer and found very nearly 
equal to 

- DS/3 . 

8 

Up to this point, the analysis has been in terms of two-dimensional angular 
distance 7 in the image plane. For comparison with experimental work reported 
in Chapter IV the desired quantity is one-dimensional variance, given by 
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where 


= <(r2 COS^ a'> _ 


r = T X + j y 

= i r cos a + j r sin a 


and 

a = reference angle . 

Isotropic behavior for optical phase tilts in the aperture plane guarantees iso- 
tropic image motion, r is independent of a , and 


cos^ Kcos^ 

= <r^> 1/2, 

The end result for <x^> is 


<x^)>-0.56 sec ^ f C^(h) dh, (58) 

Jo 


where 

D > v/XZ 

and X is expressed in radians. 

Parameter dependence is the same as for the interferometer case but the constant 
multiplier is reduced. The important results are that <x^> is independent of 
optical wavelength and only weakly dependent on telescope aperture diameter. 
Irradiance fluctuations discussed previously exhibit a rather strong wavelength 
and aperture dependence. Zenith angle and turbulence dependence expressed by 


sec d f C2(h) dh, 

where 

d = stellar zenith angle, 

predict that <x^> is proportional to air mass (sec 6) between source and observer 
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with the strength of turbulence at each altitude weighted equally. This one-dimen- 
sional result will be denoted o-^, variance of image motion, so that 

m 

For aperture diameter less than /xT the phase structure function (Tatarski) re- 
tains the same form but is reduced by a factor of about two. For this case, how- 
ever, the assumption of a constant irradiance across the aperture could not be 
used. It is likely that the resultant aperture dependence would differ from 
For D < , the inner scale of turbulence, even the phase structure function de- 

parts from ^ dependence. In the other limit as D grows very large, the 
upper bound on the phase structure function and thus <x^> is determined by the 
outer scale of turbulence. This limit is ill-defined and may be encountered as 
close as a meter or as distant as a kilometer. 

Time behavior of image motion is described statistically by autocorrelation 
analysis of a certain length sample record of data. As a first approximation con- 
sider the interferometer case discussed previously in the determination of image 
motion variance. Normalized phase difference autocorrelation, R 0 (t), for the 

d 

spatial separation, d, and variable time lag, r, is 




( 59 ) 


The denominator is by definition phase structure function for separation, d. The 
numerator may, as outlined by Kolchinskii^'^, be expanded into a sum of structure 
function terms. 
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Thus, 


~ ^B ~ ^B ^ A ^ 

\^t ^t+r ^t *^t+T ®t ^t + T ^t ^t + r/ 


(60) 


\ +^ ( 4 >A ~^b ^X 


2 '^\~^^A^T 



In order to evaluate the structure functions, certain assumptions must be made. 


(1) As in the analysis of irradiance time “dependence, consider ^^frozen-in" 
turbulence. The turbulence inhomogeneties are transferred by a wind of speed 
parallel to the line connecting point apertures A and B. 

(2) Assume a single -layer model. The turbulence is located at a distance z 
with strength C ^ and results in a phase structure function 

(d) = 2.91 k2 d5/3 2 q 2 , 

(3) The turbulence process is assumed to be homogeneous and isotropic as 
before. 

The value of phase at point A and the time A + 'T' , (t + r), is the same as that 
at the point A - Vj^ r and the time t,cj6 „ (t). Likewise, 

A- r 

^a(^^ + +r), 

0B (t + T) = ^ (t), 

and 

<^8 Ct) - + ^ (t + r) . 
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These relations allow the structure function terms in equation (60) to be expressed in 
terms of and an argument that includes r along with the point aperture 
separation, d, as follows: 


and 


or 


<;\ =D*(ldtV„r|), 


K <r) =- 


(ld+V„T|) + iD^( i d-V„ h -D*(V„r) 




E (r) = 


V r\ ^ '3 


V r \ 
N 


V r 3/3 
N 


The result is independent of strength of turbulence, distance to the turbulence, 
and optical wavelength and depends only on aperture separation, d, and the wind 
speed, Vj^. The result, of course, holds only for r values such that the phase 
structure function argument is within the inertial subrange. 

The total effect of many such layers can be examined by adding the auto- 
correlation for each layer weighted by the strength of turbulence and wind velocity 
for that layer. Thus, 
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As with the single -layer result the phase difference autocorrelation is identical 
to image motion autocorrelation. This is evident since 

0B " ^ ^ ^ ^ 

and 

^0, "Ra (r)=R^(r), 

where 

a = wavefront tilt angle 

and 

X = diffraction pattern angular shift in one-dimension. 


The Fourier transform of (''")> the angular shift autocorrelation, may be 
expressed in a form similar to equation (38) for irradiance spectral density. 

The first step is to express R^ ('^) in terms of image motion variance, . From 
equation (58) 

^00 

0-2 = .560 I (z) d z, 

m I N ' ^ ’ 

where 

D = aperture diameter ^ Z . 

Noting that 

.56 D, (R) 

0-2 

™ k2 D2 

yields the new form for R^(t) , as follows: 


= Rc 




= 1 


00 ( Id + V^rl ) + i ( Id - V^rl ) - D^(V^r) 


R<i(d) 


dz , 


or 


R ('?■) = 


. 57 

D2 a2 




_ 

cj(z) 

^0 

1 

1 

1 D + r 1 5/3 + 1 1 D _ r 1 5/3 _ 1 Vp, r 1 5/3 


(61) 
d z . 


Here the wind velocity which depends on optical path distance, z, has been brought 
inside the integration. Once again the passage from summation to an integral 
over z has been made. The Fourier transform of R (t) , 
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^00 

(f) = 4 cos (2 77 f r) (r) d r , 
can now be expressed as follows: 

w f C2 (z) F {f, V (Z), D} d 2, 

where the r and f integrations have been interchanged. This is the result 
used in Chapter V to construct turbulence profiles. 

Image Size 

Specification of image size or profile is the final step in a description 
of the atmospherically distorted stellar image. Size determination and the whole 
investigation of image quality are handled best through the use of optical transfer 
function techniques. The transfer function of an optical system, OTF, is defined 
as the two-dimensional spatial Fourier transform of the image spread function. 

In the situation examined here, optical system is composed of telescope optics 
and atmosphere. Magnitude of OTF is a measure of reduction in contrast suffered 
by each Fourier component of the object after transmission through the entire 
imaging system. It is a function of the transform variable, image plane spatial 
frequency, f. The variable, f, has dimensions of cycles per unit length. Multipli- 
cation by system focal length, D, allows spatial frequency to be expressed in 
cycles per radian field of view. This will be denoted, This, of course, is 
easily related to cycles per arc second. Relationship between OTF and image 
intensity distribution are derived in Appendix A for the particular case of star- 
light, i.e. , point source at infinity and incoherent light. The result, if T( f) 
is the optical transfer function, OTF, is 
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(63) 




T ( f) = I I (x) exp { - 2 rr L f*x)dx, 


where 

and 


f = spatial frequency vector, 
x = position vector, 

I (x) = image intensity. 


In general, T(T) is complex and can be written in terms of a modulus and 
phase as 

T ( f ) = M ( f ) exp i c;6 ( f) . 

For real, symmetric images T( f) is real and positive , and 


T ( f) . M (?) . (65) 

Analysis may proceed in terms of the modulation transfer function, M(f ) or 
MTF. Both phase and amplitude effects in the atmospherically distorted wave- 
front contribute to MTF. Phase effects can be thought of as wavefront tilts or 
bends. These lead to image tilt and blurring respectively, as previously noted. 
Amplitude effects create a random apodization of the aperture, i.e. , they weigh 
the phase effects non-uniformly. 


Just as phase and amplitude of the optical wavefront are functions of time 
so is MTF. Thus, 


M(f) - M(f, t). (66) 

This will continue to be denoted M( f) with the time dependence understood. Over 
a long period of time, as in long exposure photography, the resultant image can 
no longer be described by the telescope diffraction pattern but will be a relatively 
large blur circle. Under this condition the symmetry assumption is not too 
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severe. Time-averaging and the random character of the distortions tend to 
average out any non-symmetries that appear. For very short-exposure times, 
all effects of image motion are frozen out, leaving a sharper well-defined image. 
While this image may approach the diffraction limit it will in general exhibit a 
residual blur. A recent analysis due to Fried^^ has shown that this blur is 
primarily a focusing error. Near the optical axis this is a symmetric aberration 
caused by a quadratic term in the description of wavefront shape. Unsymmetrical 
abberations, similar to coma, are numbered third in order of importance. This 
argues in favor of retaining the symmetry assumption above. A much more 
serious question concerns the existence of spatial invariance in the image plane 
for the short-exposure case. This is the concept of an isoplanatism patch as 
mentioned in Appendix A. If such an area, in which image spread function depends 
only on the separation of two measurement points, does not exist the mathematical 
formalism breaks down and a transfer function cannot be defined. It is not 
necessary thatthe isoplanatic region encompass the entire telescope field, just 
the general vicinity of the image. This may be quite a reasonable assumption 
although further verification is needed. 


For a cascaded optical system such as the combination of atmosphere and 
telescope, the resultant MTF is just the product of individual MTF^s. Thus, 


where 



( 67 ) 


\M ( / = observed MTF averaged over a certain time, 

M (f) = telescope optics MTF which is constant in time, 
and ‘ 

= atmospheric MTF averaged over the same time as above. 

^ ^ ime 

Telescope MTF based on the assumption of diffraction-limited optics, as set forth 
in Appendix A, is 
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A. c«j < D 5 


( 68 ) 


and 


M (oS) 


2 

IT 


COS 


\ Cd 

"d" 




where 


= 0, \ O) > D , (69) 

o;= polar coordinate spatial frequency in cycles per radian field of view 
= F I f| . 


Once again the ensemble of systems necessary for a correct averaging process 
is formed by a certain time sample of data. 


Several different theoretical approaches have resulted in the same ex- 

pression for the atmospheric part of the long-term MTF, namely, 


where 




exp 



(70) 


LT = long-term time average. 

The wave structure function, denoted D (r), is the statistical description of atmos- 

w 

pheric effects on the propagating wave. It is the sum of log-amplitude and phase 
structure functions mentioned previously, i.e. , 


D^(r) =D^(r) +D^(r). 

Tatarskiis expression for (r) is very similar to (r) of equation (56). In fact, 
beyond the log-amplitude correlation distance, , they are equal. Thus, 


D (r) ^ 1 46 k2 rS/3 


^ 00 

j (z) d z, « r « y\ z , ( 72 ) 


and 


D (r) =: 2.91 k2 rS/3 


(z) d z, T Z V)Tz . ( 73 ) 
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In both cases the assumption is made that 


/X z » , 

As Fried has pointed out, the conditions on r may be replaced by the same con- 
ditions on telescope aperture diameter, D. The last equation applies in general 
for telescopes of diameter 


D^IO centimeters, 

and in particular for the telescope used in the observations reported in Chapter 
IV, Spatial variable, r, in the aperture can be expressed in terms of image 
plane spatial frequency, oj, by 

r = k a) . 

The final desired form for D (u') is 

W ' 


(co) =2.91 (4 7T^) sec ^ f ^ ^ 

Jq 

where C is a constant for a particular set of data and expresses the parameter 
dependence. The long-term atmospheric MTF becomes 




= exp 



Spatial frequency dependence is contained in the factor. Image size de- 

pendence on C is revealed by applying the Fourier transform result of Appendix A, 
so that 


I (s) 


2 77 


2 77 



(2 77 s cj) d o; 


Jq (2 77 s co) d . 


(75) 
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Here for convenience a; is expressed in cycles per arc second and distance, S, 
in the image plane is expressed in arc seconds. The 5/3 power law is close to a 
square power law and this approximation will be made. * When this is done the 
rather complicated integral above can be easily evaluated as follows; 

p CO 

I(s)~27tJ cj exp (2 tt s oi) d o) 


2 77 

= exp 

C 


C ^2 77 


This is a Gaussian shape with variance, o-^ , such that 

Cj2 C "4 772 
s 


(76) 


= 2.91 sec 0 


^00 

C2 (h) dh . 


Here is seen to be proportional to the negative one-third power of wavelength, 
secant of zenith angle, and the integrated strength of turbulence. This is essentially 
the same result obtained for variance of image motion, a2 , with the weak 
dependence on aperture replaced here with a weak dependence on wavelength. 


Short-term MTF has also been treated in the literature. An expression that 
is similar to that discussed by Fried^^is 





1 - 1.026 



(77) 


*See Appendix C. 


= exp 



1.026 
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where 

d^VTT 

and »co 

C 2.91 (4 7/2) \-i/3 sec (9 (h) d h . 

•^0 

The factor is the same as that found in long-term MTF and dominates the be- 

havior of <^m(a;)^^g^ for a; < < ( VD)o Inclusion of the second term with its 
factor modifies the behavior for large spatial frequencies, ^ < (A./D). The factor 
\/D is the resolution cutoff spatial frequency for a telescope of diameter D. This 
may be observed by substitution in equation (68). Since short-term MTF depends 
on the particular telescope diameter used, ^ tvne atmospheric 

MTF. As Fried has pointed out, an MTF can be assigned to atmosphere only in 
the long-term case. 

Long-term MTF can be measured experimentally. The length of time must 
be long enough to include all the effects of image motion, however. This implies an 
averaging time of tens of seconds, while the typical MTF measuring device com- 
pletes a scan in a fraction of a second. The experimental device described in 
Chapter V employs instead an average measurement of many short -exposure 
images. This data can then be combined with a simultaneous measure of image 
motion, through the analysis presented below, to achieve an estimate of long- 
term MTF. This is the desired result since it is independent of all lens param- 
eters and may truly be assigned to the atmosphere. 

Long-term MTF can be defined in terms of short-term MTF and a transfer 
function assigned to image motion by 

<“ (")«>„ ' <">«>„ M (">.>„ • (78) 

If image motion and short-term blur can be considered uncorrelated, it follows 
that 
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3 7 

Following the lead of Hufnagel the image motion MTF can be expressed in 
terms of one-dimensional image motion variance, 0*2 , by 

m 



( 79 ) 


where is determined from a sample of data whose length roughly corresponds to 
long-term time average. 


A rough check on this procedure can be made by assuming that « 
can be adequately represented by equation (77) for <M (^)y^>gy alone. Here the 
averaging over the long term is applied to the wave structure function appearing 
in the exponent of <M (a;) > . Thus, 

^ ST 

<« ' <« <“ . 
or 


exp 




± D (u>) 

2 w ' ' 


1 - 1.026 


D 




(80) 


or 





. 513 




but 


Hence, 


f co 

(h) dh . 



~ exp 


4 7T^ cl)^ 1 . 49 sec 0 



(h) dh >. 


(81) 
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Contents of the exponential are within a constant factor equal to the expression 
for one-dimensional image motion variance, cr^ , derived earlier. 

Accurate determination of long-term MTF affords the observer a unique 
opportunity to specify optical performance of the atmospheric channel. Here 
with one result, independent of the optical system used to determine it, resides 
a powerful tool for predicting the effects on larger and more complicated systems. 
In the following two chapters, the experimental part of this research will be dis- 
cussed and data presented. An attempt will be made to construct a turbulence 
profile. The strength and shape of this profile remains as the largest single un- 
known in all of the proceeding theory. 
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CHAPTER IV 



ACQUISITION OF EXPERIMENTAL RESULTS 
Stellar Image Monitor Apparatus 

Experimental data for this thesis research is provided by a device developed 
by the author and fellow workers at Goddard Space Flight Center during 1968 and 
1969. It was designated the Stellar Image Monitor and was designed to provide 
quantative data for optical site selection. The Stellar Image Monitor, or SIM, 
is based on a design published by Ramsay^® and a technique described by 
Lindberg 3^. Ramsay performed some experimental work with starlight and 
later Coalman used the device extensively for horizontal optical propagation 
studies. SIM hardware produces three analog voltages containing information on 
image intensity, size, and motion. These voltages can be analyzed for the desired 
statistics of the stellar image. The heart of the device is a modulation transfer 
function measurement performed by square -wave modulating or chopping in the 
telescope focal plane. Chopping is performed at a selected but variable spatial 
frequency. Operation at various spatial frequencies, detection of the chopped light, 
and electronic processing yield various points for an MTF curve. Image intensity 
is monitored by low -pass filtering of the chopped signal. Image motion information 
is derived from electrical phase comparison of signals generated in chopping the 
stellar image and a reference light source. The Goddard device differs from the 
original Ramsay design in the way spatial frequency is varied, in electronic 
processing, and in data collection and reduction. Data in the form of analog 
voltages is recorded in the field on magnetic tape. Off-Hne processing by digital 
computer results in a statistical record of the recorded data, facilitating compari- 
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son with theoretical analysis presented in the previous chapters. Data results 
presented in this Chapter were obtained between fall 1968 and spring 1969 while 
the SIM was being field tested at the Goddard Optical Research Facility. The 
device is presently in an operational program of site evaluation at the Smithsonian 
Astrophysical Observatory at Mt. Hopkins, Arizona. In reporting results, a 
concentration will be placed on those deemed most effective in estimating strength 
of turbulence along the vertical path and in testing the theory. Figure 3 illustrates 
the experimental concept and Figures 4 through 6 are photographs of the equipment 
used. 



X STARLIGHT 
// 


TURBULENT ATMOSPHERE 




Figure 3. Stellar image monitor experiment concept 
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Figure 5- Detail of electronic data collection equipment 
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In the SIM, starlight is collected by a 6-inch diameter doublet lens system 
of 48 “inches focal length. Figure 7 indicates the optical system components. As 
converging starlight enters the experiment package, it encounters a partially 
silvered mirror. This surface, acting as a beamsplitter, allows 30 percent of the light 
to pass through and reflects 60 percent. Transmitted light enters a pentaprism 
and is directed through a corrector plate to focus on a ring and crossline reticle. 

The field-of-view corresponding to the ring^s diameter is 14 minutes of arc. 

Observation of the stellar image on the reticle allows pointing of the telescope 
to insure on-axis operation. Light reflected by the beamsplitter is reflected 
again by a front surface mirror and allowed to come to focus. A microscope 
objective serves as relay optics to focus an enlarged stellar image on a rotating 
glass disk. An eyepiece and mirror assembly mounted in a sliding tube are used 
to visually check image quality and insure proper alignment of the image on the 
glass disk. The disk is composed of alternate clear and opaque pie-shaped 
sectors as detailed in Figure 7. Each sector subtends one-half degree. The 
pattern was produced at Goddard by photographic -reduction and etching from a 
machine-drawn original. The disk is belt driven from a D.C. motor-servo system 
at a constant rate of 517 rpm. Speed regulation is better than 0.5 percent. Light passing 
through the disk is thus square-wave modulated with a time frequency of 3.1 KHz, 
determined by motor speed and number of sectors on the disk. The spatial fre- 
quency of chopping is inversely proportional to the distance from the center of the 
disk and is determined by sector width at a given distance. Since one clear and one 
opaque sector determine one cycle of chopping, the spatial frequency in cycles 
per millimeter is the reciprocal of the width of two sectors. Spatial frequency 
can be made independent of lens parameters by multiplication by system focal 
length. Cycles per millimeter become cycles per radian field-of-view or cycles per 
arc second. Starlight is brought to a focus one inch from the center of the disk. 
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Figure 7. Optical system and detector package 


For unity power relay optics the spatial frequency of chopping is 2.26 cycles per 
millimeter or 0.0134 cycles per arc second. Higher powers effectively increase the 
spatial frequency by enlarging the image on the disk. Below are listed the linear 
increases in spatial frequency achieved with various powers of microscope 
objective. 


Light transmitted by the disk passes through an optical bandpass filter and is 
received by a photomultiplier tube. The photomultiplier tube is an EMI Model 
9558B with S-20 response. The filter is a Corning glass filter, Number CS 3-71. 

o o 

The combined spectral half width is approximately 1200A centered at 5000A 
(Figure 8). 
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RELATIVE RESPONSE 


Table 1 


Microscope Objective Power and Spatial Frequency 


Microscope Objective Power 

Spatial Frequency 
(Cycles/Arc Second) 

1.0 

0.0134 

2.6 

0 . 0348 

3.5 

0.0469 

4.0 

0.0536 

6.0 

0.0803 

10.0 

0.134 

21.0 

0.281 

45.0 

0.602 

60.0 

0.803 

98.0 

1.31 



Figure 8. Combined relative spectral response 




Light from a small D.C. powered lightbulb is also focussed on the disk to 
act as a phase reference source. A thin glass plate in the converging light from 
the bulb can be tilted to adjust the position of this image on the disk. Chopped 
light is received by a PIN-10 photodiode. A field effect transistor preamplifier 
and biasing network prepares the photodiode signal for entry into the electronic 
preprocessing system. The photomultiplier signal is developed across a one 
megohm load resistance and unity gain preamplifier stage. 

Electrical signals from the two detectors are fed into an electronic pre- 
processing system (Figure 9). The main function of this device is to produce 
three voltages containing information on image intensity, modulation index, and 
motion. Operational amplifiers are the primary circuit constituents. The 
entire instrument is D.C. powered by a regulated power supply providing 27 
volts at approximately 0.2 amps. Solid state voltage regulators provide the 
plus and minus 12 volts required for operational amplifier biasing. In addition, 
a programmable D.C. -to-D.C. converter provides up to 2000 volts for the 
photomultiplier tube. 

The positive going input from the photomultiplier preamplifier stage is in 
the form of an amplitude -and frequency -modulated square wave at the chopping 
frequency, 3.1 KHz. This signal and the photodiode signal are not perfect square 
waves since image size is not negligible compared with sector width. Only the 
fundamental sine wave component of the square wave is of interest, however. 

As shown in Figure 9, the photomultiplier signal proceeds through two inverting 
D.C. amplifier stages. The first has a gain of 10 and the second incorporates a 
variable gain from one of 20. At this point, the amplified signal proceeds 
two ways. The first path is through a low -pass filter of the multiple feedback 
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SYSTEM DETAIL 

Figure 9. Electronic preprocessing system 












design having a corner frequency of 200 Hz. The filter's output is just the meas- 
urement of total image intensity and its scintillation. This D.C. or low-frequency 
information output 1, is also used in the determination of modulation index. The 
second path is through a multiple feedback bandpass filter centered at the chop- 
ping frequency. This frequency belongs to the fundamental sinusoidal component 
in the square wave. The 3-db bandwidth of the filter is plus or minus 200 Hz 
around the center frequency. An analog divider receives the outputs of both 
filters and performs the operation A.C. divided by D.C. The output of the 
divider is the normalized chopping frequency component, that is, the strength 
of the component with respect to D.C, As the image size grows larger due 
to turbulence or as spatial frequency is increased, the A.C. component 
suffers a reduction. Size changes affect only the A.C. component while inten- 
sity changes or scintillation are present in both A.C. and D.C. components 
and cancel out in the division process. The divider's output is received by a 
peak detector which produces a D.C. voltage, output 2, proportional to the 
modulation index or relative strength of A.C. to D.C. Fluctuations in modula- 
tion index can be followed from D.C. to the corner frequency of 200 Hz. If the 
peak value of A.C. is equal to the D.C. signal strength, the modulation index is 1.0. 

Motion of image center of gravity is determined by phase comparison of 
signals at the chopping frequency, 3.1 KHz, from the two detectors. As the stellar 
image moves on the chopping disk, it produces a waveform shifted in phase with 
respect to that from the phase reference source with its stationary image. The 
preamplified reference signal from the photodiode is first sent through a bandpass 
filter identical to that used in processing the phototube signal. The outputs of 
both bandpass filters are sent through limiter stages to provide good, clean 
square waves. An integrated-circuit phase comparator receives the square 
waves and outputs a pulse train whose duty cycle is proportional to the 
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relative phase of the signals. An RC averaging network with response from 
D.C. to the corner frequency of 200 Hz measures the duty cycle. Relative phase 
of the two signals is just the relative position of the two images on the disk along 
a direction perpendicular to a disk radius. Thus, output 3 becomes a voltage 
proportional to one-dimensional motion of stellar image center of gravity. 

The results of electronic and optical calibration tests of the preprocessing 
system are given in a following section. These tests were necessary to insure 
proper operation of the electronics and to relate the voltage outputs to parameters 
specifying image quality. 

In order to preserve the output data of the electronic preprocessing system 
in a form suitable for data analysis, recordings are made on magnetic tape. For 
this purpose, an Ampex FR-1300 portable instrumentation tape recorder was 
incorporated into the system. This is a seven track machine equipped in a 
standard TRIG configuration for one-half inch wide magnetic tape. The three outputs 
of the electronic preprocessing system have bandwidths from D.C. to 200 Hz. 

The D.C. requirement means a frequency -modulation (FM) recording scheme must 
be used. Time signals are also recorded in order to properly characterize and 
identify the data. The NASA one-per-second and one -per -minute binary time codes 
are the ones used. These can be recorded in the direct (AM) mode. Recording is 
done at a tape speed of 7.5 inches per second. The FM and AM record bandwidths 
are from D.C. to 2.5 KHz and from 50 Hz to 38 KHz, respectively. A 10 KHz 
sine wave is also recorded in the direct mode as a reference signal. Tape speed 
errors made in data taking can be eliminated during automatic data processing 
by using the 10 KHz reference signal. Such errors would result from fluctuations 
in power-line frequency. Signal-to-noise ratio specifications for the tape recorder 
are greater than 40 db at the tape speed used. 
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Operation of the Stellar Image Monitor Experiment involves nighttime stellar 
observations at various zenith angles, azimuths, and different times of night. 

Two operating personnel are required, one to adjust telescope pointing and op- 
tical alignment and the other to monitor the electronics and operate the tape recorder. 
In acquiring a star, successive use is made of a wide field viewfinder and narrow 
field, high -power eyepiece. This insures operation of the telescope on-axis. The 
focussing and centering microscope objective holder is then used to bring the 
stellar image to the correct focus and position on the chopper wheel. This process 
is monitored by the sliding eyepiece. After proper alignment is achieved, the 
eyepiece is removed and starlight enters the detector package. Once the photo- 
multiplier high voltage control and D.C. amplifier gain have been adjusted for the 
proper intensity output voltage, the instrument is ready for data recording. 

Data -taking procedures are tailored closely to the demands of data reduction. 

To this end, calibration voltages, a star identification code, and background runs 
are recorded on magnetic tape along with the data, time, and reference signals 
(see Figure 9). A voice-edge track is used to note experimental parameters such 
as spatial frequency, telescope aperture, photomultiplier tube high voltage, and 
amplifier gain. Data taking is divided for convenience into four distinct categories. 

A tape is the largest category. It is the contents of one reel (2,500 feet) of mag- 
netic tape and may include data from several nights , the second category. The 
third level is called a star . It includes all data taken in succession with a particu- 
lar stellar source. Included here may be from one to ten data runs , the final 
category. 

A data run begins when the tape recorder is turned on and ends with tape re- 
corder off. It is characterized by a three-digit data run number and a particular 
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spatial frequency. Nominal operation would include five data runs per star, 
each conducted at a different spatial frequency. Below is listed the format of 
a data run. An example is given in Figure 10. 


Table 2 

Data Run Format 


Event Number 

Event 

Length (Seconds) 


START RECORDER 


1. 

Tape leader at ground potential 

20 

2. 

Positive 6 volts calibration 

10 

3. 

Ground potential 

5 

4. 

Negative 6 volts calibration 

5 

5. 

star number himdreds digit 

3 

6. 

Star number tens digit 

3 

7. 

Star number units digit 

3 

8. 

Data 

10 to 600 

9. 

Defocus data 

5 

10. 

Background data 

5 


STOP RECORDER 



The first seven events are fed automatically to the tape recorder by a com- 
mutator arrangement. The commutator selects the desired D.C. voltage level 
for the specified length of time. Tape leader allows time for the tape recorder 
transport mechanism to come up to speed and stabilize. Known calibration voltages 
are recorded to provide a voltage standard in data reduction. The three star 
numbers are the hundreds, tens, and units digits that comprise the number as- 
signed to a particular star. This number may be found in Apparent Places of 
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Z/ZZ 


DATA 


STAR 

TAPE LEADER CALIBRATION NUMBER 


DEFOCUS BACKGROUND 




TIME SCALE (SECONDS) 


10-^5*f5*|3 


10 TO 600 


5+H 


Figure 10. Data run format 



Fundamental Stars . Each digit is a certain percentage of plus calibration volt- 


age. For example, the star number 297 would comprise voltage levels 20, 90, 
and 70 percent, respectively, of plus calibration. Actual star data is taken from 10 
seconds to 10 minutes. Defocus data is generated when the telescope operator 
dehberately defocuses the stellar image on the chopper wheel. This is used to 
characterize the modulation-index data. Background data requires the operator 
to point the telescope so that the star is just beyond the field-of-view. Thus, all 
the effects of stray light, sky background, and electronic offset can be recorded. 

Data from the three outputs of the SIM is recorded simultaneously in each 
data run. Time and reference signals are recorded continuously without the 
calibration voltages and star code. Below is a listing of the seven tracks of the 
magnetic tape recorder, their functions, and modes of recording. 


Table 3 

Tape Recorder Format 


Track Number 

Function 

Record Mode 

1 

Spare 

FM 

2 

10 KHz reference 

AM 

3 

Image motion data 

FM 

4 

! 

Intensity data 

FM 

5 

Modulation index data 

FM 

6 

BCD time code 36 bit 

AM 

7 

BCD time code 28 bit 

AM 


The electronic equipment described above is mounted in a standard 5-foot 
instrumentation rack as shown in detail in Figures 5 and 6. The following items 
are included: 
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(1) Electronic preprocessing'systein . 

(2) Calibration voltage commutator. 

(3) Tektronix type 422 oscilloscope. 

(4) Newport Labs series 200 digital panel meter. 

(5) Trygon HR40-3B D.C. power supply. 

(6) Electro-Craft E-650 motor speed control system. 

(7) Ampex FR-1300 instrumentation tape recorder. 

The oscilloscope and digital panel meter are used to monitor signal levels during 
data recording. Time and reference signals must be fed into the rack from an 
external time code generator. 

The three outputs of SIM are randomly fluctuating voltage levels representing 
optical fluctuations in the stellar image being monitored. Similar to the atmos- 
pheric mechanisms responsible, these optical fluctuations have a large dynamic 
range and periods as short as several milliseconds. In Figure 11, sample chart 
recordings of each of the three data tracks are presented. Clearly, one must 
turn to statistics in order to characterize the data. Consideration must be given 
to averages, variances, correlations, and power spectrums. 

Once data has been recorded on magnetic tape, the first step in data reduc- 
tion converts the information into a digital format. For this, the facilities of the 
Information Processing Division at GSFC are used. Analog data tapes are played 
back through an analog to digital (A/D) conversion unit and digital tapes are pro- 
duced. A/D sampling is at a rate of 4 KHz. A sample is taken every quarter 
millisecond. The data tracks are sampled sequentially with the net result that 
all of the three outputs of the SIM plus a spare fourth track are sampled every 
millisecond. Fundamental sampling theorems require that the sampling rate be 
twice the highest frequency component in the data. This insures correct 
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representation in digital form of the original data. Sampling 1000 times per second 
allows the data to contain frequencies as high as 5000 Hz. This is more than 
enough to contain the 200 Hz bandwidths of the SIM electronics and, from pre- 
vious experience, the expected bandwidths of optical fluctuations. The 36-bit 
one -per -second time code is also digitized and used primarily as a bookkeeping 
device. The digital tape contains information in binary coded decimal format and 
is ready for computer processing. 

Processing is done on an IBM system 360/91 computer through programs 
written in the Fortran IV language. The following description of data analysis 
concentrates not on specific program instructions but rather on the method of 
statistical analysis. This method conforms to standard practices in statistical 
analysis. The method is best examined by following a sample data run. 
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The first information input to the computer for each data track is a series 
of samples of the positive, zero, and negative calibration voltages , v^ , and 
v^. As mentioned previously, v^ = +6.0 volts, v^ =0.0 volts, and v^ = -6.0 volts. 
The computer is programmed to determine average values , X 2 , and x ^ of 
the samples for each calibration voltage. A linear relation is assumed to exist 
between actual voltage v and sensed value x, that is, 

V = A (x - B). 

Two of the resulting three simultaneous equations can be solved for scale 
(A) and tare (B) factors. This calibration is used later in the data run to assign 
true voltage readings to the data. It serves to eliminate inevitable gain changes 
and D.C. offsets present in the recording and reproducing process from magnetic 
tape. The computer is next presented with samples of the three digit star num- 
ber. Once again, average values are determined and are expressed in percent- 
ages of positive calibration. This serves to identify the star number. The 
computer has stored a list of over 100 stars according to star number, 
right ascension, and declination. Use of the following information allows calcu- 
lation of stellar zenith angle and azimuth: 

(1) Time of data run. 

(2) Sidereal time at Greenwich, England, for zero hours Universal Time on 
the day of observation. 

(3) Stellar coordinates, right ascension and declination. 

(4) Longitude and latitude of observation site. 

Data samples then begin to enter the computer. If 20 seconds of data is to 
be analyzed, the computer will be presented with a series of 20,000 data points 
for each of the three tracks. All statistics subsequently calculated from the data 
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points will be time averaged over the length of the data run (20 seconds). It is 
necessary to average data over a time long enough to assure good statistical 
accuracy. Estimates based on equivalent bandwidth and resulting number of in- 
dependent samples suggest at least a 10 -second averaging time for SIM data. 
On the other hand, if averaging time is too long, atmospheric or other experi- 
mental conditions can change significantly, destroying the meaning of the data. 
Twenty seconds seems to be a reasonable compromise. 


Each data point is assigned a voltage level, v, according to scale and tare 
factors determined by calibration. Method of analysis is then split into two main 
procedures. In the first procedure, the voltage range that includes all the data is 
divided into a large number, K, of class intervals. The width, Av, of the class 
interval, i, may be chosen at will but has as its effective lower limit the voltage 
resolution of A/D equipment. The following are t3T3ical parameter values: 


and 


voltage range = 0.0 to 10.0 volts, 
K= 500, 

Av = 0.02 volts. 


The computer tabulates the number, n(v. ), of occurrences of the data record in 
each class interval and also the percentage occurrences, p(v. ). This in effect 
establishes a probability density histogram, P(v. ), representii^ in digital form the 
probability density function, P(v) , by 


P (Vj^) 


n (v.) 


where N = total number of data samples. Then, 


P(v.) 


P (v.) 
A v 
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and 


K 

2^P(v.)Av = 1 . 

i=0 

The cumulative density histogram is related to P(V| ) by 


i=l 


This represents the total probability for voltage less than or equal to the given 
voltage, Vj . Probability density and cumulative density are related theoretically 
by 


and 


P(v) 


dCCv) 

dv 


C(v) 



P( V) dv . 


The computed probability density function is used to determine means, variances, 
and higher moments of the data record. Thus, 

N 


Mean - v 




Av 


i=l 


j moment about the mean 


IN 

2^(v. - v)jp(v.) Av 
i=l 


m-* 


where j = 2 is the variance, The sample mean and sample moments thus com- 
puted are not unique for a given sequence of N data points, since they are determined 
by choice of voltage interval endpoints and by number of class intervals. 
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The parameters coefficient of kurtosis, k, and skewness, s, help to charac- 


terize the probability density function. These are given by 

.4 


and 


k = 


s 


m" 




- 3.0 




(m2) 


3/2 


Coefficient of kurtosis relates to the width of the distribution and skewness re- 
lates to its symmetry. For a Gaussian distribution, both k and s will be zero. 
Random voltage statistics can best be summarized in one parameter, coeffic- 
ient of variance, where 

Coefficient of variance - =i • 

V 

This parameter has the form of a percentage modulation. It is the ratio of stand- 
ard deviation of the fluctuating signal to its average value. 

The second main direction of analysis is based on data time series, the suc- 
cession of data samples maintained in their original positions as a function of 
time. The basic tool in this analysis is the sample autocorrelation function, 


1 /-T-r 

R(t) - j ~ I v(t) v(t + r) dt , 

Jo 

where 

r = lag time , 

T = length of data sample, 

and 

v(t) = voltage level at time t. 

The autocorrelation function for random data describes the general dependence 
of the values of the data at one time on the values at another time. An estimate 
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of the autocorrelation between values of v(t) and times t and t + r is just the 
product of the two values averaged over the length of data sample. Here, r is the lag 
time or separation time. For computer analysis, the discrete representation is 
used. Thus, 


where 


N- r 


R(r) 


- RC'CAt) - v/ 

i = 0 


'i = lag number = 0, 1, 2, ... m, 


and 


V 


= V. 
1 


m = maximum lag number, 

N = total number of samples, 
t = interval between samples, 

“ V = voltage level for the i^*^ sample minus the average voltage level. 


For a good estimate of R(t) m should be chosen one tenth of N or less. Typical 
parameters are 


and 


N = 20,000, 
m = 200, 

At = 1 millisecond , 


Although R(r ) is a useful quantity by itself, it is more useful when applied to an 
estimate of spectral density, G(f), where 


G(f) 


1 im 
Af-"0 


(f) 

Af 


Spectral density can be thought of as voltage variance per unit bandwidth. In order 
to proceed from the time domain, R(r), to the frequency domain, G(f), we must use 
the Fourier transform of R(r)with respect to r, so that 
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G(f) 


R(r) cos (2 tt fr) dr . 



In discrete form for computer analysis, this becomes 


G(f) 


m 

4At R(^At) 
^=0 




COS 


where 

f = frequency (Hz) 

and 

= 500 Hz = frequency cutoff of data. 


The equivalent bandwidth or frequency resolution in this analysis is , where 


B 


e 



For the previously assumed parameters, B^= 5 Hz. The raw estimate spectral 
density requires smoothing because the variability of the estimate does not de- 
crease with increased sample size or record length. Smoothing or weighting 
nonuniformly can be accomplished by multiplication by the Hanning lag window, 
D('tAt) , 


Thus, 


D(^At) 




COS 


'TT'l 

m 


G(f) 



D(^At) R(^At) cos 



For graphical presentation of spectral density, normalization to unit vari- 
ance may be employed. The variance of the signal in each frequency interval is 
divided by variance for the full signal. This insures that all spectral densities 
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will have the same area under the curve regardless of the strength of turbulence 
fluctuations. This method aids in plotting strong and weekly fluctuating signals 
on the same graph. The unnormalized plot should be used to facilitate compari- 
son of strength of signal in a given bandwidth for many different data runs. 

Although both methods of statistical analysis apply in general to all data 
from the SIM, each data track must be handled in a particular way. The image 
intensity track represents irradiance fluctuations of starlight. Voltage produced 
by the monitor is directly proportional to irradiance, I, optical power per unit 
area at the receiver. Before the connection can be made, the effect of optical 
and electronic background must be removed from voltage statistics. The pri- 
mary effect of background is the introduction of an unwanted D.C. offset into the 
data. Computer measurement of background data taken at the end of each data 
run allows the D.C. offset to be subtracted before calculations are made. The 
dimensionless parameter coefficient of variance, here called CIV, coefficient 
of irradiance variance, has been the most widely used measure of irradiance 
fluctuations. The standard format is a plot of log^ [ CIV ] versus log^ [secant 
(zenith angle)] of the stellar source. This plot is characterized by the diameter 
of the telescope aperture and, of course, the turbulence state of the atmosphere. 

As indicated previously, log-amplitude and, thus, log-irradiance behave 
in a Gaussian manner. The logarithmic behavior of irradiance is evidenced by a 
low average value compared to high peaks, often many times the average value but 
of short duration. Thus, image intensity could be discussed in terms of log- 
amplitude statistics. Log- amplitude can be written as one-half the natural log 
of the ratio of instantaneous irradiance, that is, one-half log-irradiance, so that 
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l/i is just equal to v. /v since the proportionality constant cancels. Thus, 

1 

= "o' log • 

Statistics of log amplitude can be computed in the same manner as voltage sta- 
tistics with the substitution for v. . Log amplitude variance, , now becomes 
the central parameter. Summary information is a plot of log amplitude variance 
versus zenith angle. Spectral densities of both irradiance and log amplitude are 
particularly important for they reveal much about the atmosphere through their 
demonstrated dependence on zenith angle and upper altitude wind velocity. 

The modulation-index track results in a voltage proportional to the modulus 
of the optical transfer function for the particular spatial frequency used in data 
taking. Voltage statistics are appropriate and require only a conversion factor 
to modulus statistics. Laboratory tests generate a plot of modulus versus modu- 
lation index voltage. Modulus is just the depth of modulation of the square wave 
input to simulate the waveform from the detection of chopped starlight. One- 
hundred-percent depth of modulation corresponds to a modulus of 1.0. As the 
depth of modulation is varied, the analog divider and peak detector respond with a 
change in voltage. Modulus values obtained from this curve are corrected for 
the effect of electronic and optical background. This is done with average 
voltage and average background data obtained from the image intensity track. 

The average value of defocus voltage serves as the zero reference for modula- 
tion index voltage. The variance of modulus is read directly from the calibration 
curve. Data reduction is incomplete without the time series information provided 
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by autocorrelation function and spectral density. This analysis can be performed 
directly on modulation index voltage. Secant of stellar zenith angle, azimuth, 
and spatial frequency serve to characterize the data run. 

The image motion track may be handled entirely in terms of voltage statistics, 
since optical phase, image motion, and output voltage of the electrical phase com- 
parator enjoy a linear relationship. Image motion in the telescope focal plane 
is the result of random deviations or tilts in the angle of arrival of various parts, 
or all, of the incoming wavefront. These tilts are actually optical phase differences 
following a Gaussian distribution. Stellar image motion with respect to reference light 
image causes a change in the relative electrical phase of waveforms from the 
two detectors. Suitable laboratory calibration yields a conversion factor in volts 
per arc second for each spatial frequency. Voltage statistics are thus converted 
to image motion statistics. Present in the image -motion voltage output is a linear 
drift attributable to telescope tracking error. This is compensated by a computer 
generated least squares fit of a straight line to the data. The conversion to 
average image position is meaningless in data reduction, since it is not important 
to know precisely where the image is but how it fluctuates with respect to a fixed 
point. The important parameter is image -motion variance, . Either this 
number or rms value, cr^ , is plotted versus secant zenith angle of the star. The 
measured motion is only the component perpendicular to the chopper disk radius. 
Spectral densities of image motion can be scaled directly from voltage spectral 
densities. 

At the conclusion of computer analysis of stellar data for a particular night, 
a data summary sheet is prepared. The following information is tabulated for 
each data run. In addition, data for different stars is combined to yield plots of 
logg (CIV) versus log^ [secant (zenith angle)] and phase variance versus secant 
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(zenith angle). A modulation transfer function is plotted for each star observed 
at four or more spatial frequencies. This places the data in a form suitable for 
interpretation. 

Table 4 

General Parameters for Each Data Run 


1. 

Data run number 

2. 

Star name and number 

3. 

Zenith angle 

4. 

Secant (zenith angle) 

5. 

Logg [secant (zenith angle)] 

6. 

Azimuth 

7. 

Time interval for data run 

8. 

Spatial frequency 

9. 

Aperture diameter 


Table 5 

Statistical Summary for Each Data Run 
INTENSITY TRACK 

1. CIV = V = coefficient of irradiance variance 

2. cr^ = log amplitude variance 

3. V = signal mean 

4. log e (CIV) 

5. Vg = background mean 

MODULATION INDEX TRACK 

2. V = signal mean 

3. Vg = defocus mean 

IMAGE MOTION TRACK 

1. V = signal mean 

2. = motion variance 

3. - motion rms 

m 
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Prior to deployment of the Stellar Image Monitor Experiment in the field, 
a series of laboratory tests and calibration runs was performed. The first 
system component to be tested was the electronic preprocessing system. 

Figures 12 and 13 show measured frequency and phase response of the bandpass 
and lowpass filter sections. Bandwidths are 3100 ± 200 Hz and 200 Hz, respectively- 
Frequency response was measured by using a sine wave signal generator in place 
of inputs from photomultiplier and photodiode. Phase response was measured on 
an oscilloscope. Filter input and output were applied to the vertical deflection 
and horizontal sweep respectively. The resulting lissajous pattern displayed 
on a phase reticle gave readings of phase difference. Gain and linearity tests 
of the amplitudes were also performed. This served to calibrate the image 
intensity track. 



83 




Modulation index calibration required a measure of combined response of 
analog divider and peak detector. Sine waves at the chopping frequency were 
input in place of signals from the photomultiplier. Modulation index (peak A.C. 
divided by D.C.) of this input was measured at the input terminals of the analog 
divider. Measured values of peak detector output voltage for the given modulation 
index or percent modulation are shown in Figure 14. This curve applies for D.C. 
voltage in the range -1.0 to -3.0 volts. It should be noted that linearity is not an 
important consideration in this output since the expected fluctuations will be 
small compared to their average value. Percent modulation determined from 
this curve must be multiplied by 


D.C. + BOND / d.C. GAIn \ 
D.C. \A.C. GAIN^* 
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Figure 14. Laboratory calibration of modulation index 


to correct for filter gains and any D.C. offset voltage. The symbol BGND, or 
background, represents this offset. D.C. + BGND information is available from 
the image intensity track. Filter gains were measured during calibration and 
their ratio is 


D.C. GAIN 

A.C. GAIN ■ 1 025. 

Calibration and testing of the image motion electronics employed an artificial 
star source, Figure 15. A point source at the focus of an f/10, 20 -inch diam- 
eter collimating mirror represented the small angular size and parallel light of 
a star. The point source was either white light or 0.6328 micron red light from 
a Spectra Physics 131 helium-neon laser. Light from an incandescent lamp 
filament focused on a variable size pinhole served as the white light point source. 
In both cases a further reduction in spot size was achieved through microscope 
objectives. 35 



Each thousandth of an inch in the focal plane of the 20 -inch collimator 
corresponds to one arc second. Source sizes were kept smaller than the theo- 
retical angular resolution of 6 -inch diameter optics. 

Collimated light produced by this test setup was allowed to enter the 6- 
inch diameter optics of the Stellar Image Monitor after traversing about 40 feet 
of laboratory air space. The point sources were mounted on an optical bench 
equipped with a translation slide. Movement of the point source in increments 
of 0.001 inch in the focal plane could be measured there and produced the cor- 
responding number of arc seconds deviation. The D.C. voltage output of the SIM 
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phase comparator was recorded as a function of this arc second displacement. 
The output of the phase comparator was not a pure D,C. voltage but contained 
frequency components due to room vibration, air path turbulence, and chopper 
pattern centering error. The resultant plot. Figure 16, yields the desired cali- 
bration factor in volts per arc second. This factor depends linearly on spatial 
frequency. For the spatial frequency used in the calibration the rms error was 
approximately 0.1 arc second. It should be noted that the calibration curve is 
cyclic with a linear response in the mid-range of voltages. Thus in actual data- 
recording the average value of phase comparator output should be maintained at 
approximately 3.5 volts. A check on spatial frequency was made by noting the 
number of cycles of the response curve per arc second. The results agreed with 
microscopic objective power being used and calculation of chopper wheel sector 
width. 



Figure 16 . Laboratory calibration of image motion 
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Using the artificial star system described above optical transfer functions 
for the SIM in white light and red light were generated (Figure 17). This was 
done by varying microscope objective power in the SIM. Measured OTF shows 
a degradation of performance in white light while monochromatic red light per- 
formance is close to the diffraction limit. Chromatic aberration causing this 
effect must be present in the SIM optics since the collimating system is reflective. 



SPATIAL FREQUENCY (CYCLES /ARC SECOND) 
Figure 17. Optical transfer function modulus 


Error Analysis 

In examination of possible error sources to account for inaccuracies in the 
analyzed data and spread in data point values, three major factors have been de- 
termined. Listed in the order of importance of their effect on the data, they are — 

(1) Nonstationarity of the atmospheric processes. 

(2) Statistical sample size. 

(3) Electronic signal -to-noise ratios. 
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To even the most casual observer of optical phenomena resulting from 
atmospheric turbulence, the most striking fact is their nonstationarity . Not 
only do these phenomena vary from night to day and season to season but over 
periods as short as several minutes, large departures from stationarity may be 
observed. In order to characterize the phenomena, the scientist must take 
samples of finite length and perform statistical analyses. He is limited on the 
lower end of the time scale by his desire to get enough data to form a reasonable 
statistical sample and on the upper end by nonstationarity effects and the fact 
that his equipment can only handle so much data. Estimates of the magnitude 
of nonstationarity are in dispute but it is generally conceded that the problem is 
much worse in the ground boundary layer of the atmosphere. Here local terrain 
effects become dominant. Lawrence, Ochs, and Clifford have reported fluctua- 
tions in the ground strength of turbulence, C^, by factors of 10 to 100 
over periods of a few minutes. Image size and image motion data are mostly de- 
termined in the lower layers and thus may be expected to exhibit the greatest 
nonstationarity. Irradiance, determined by higher altitude effects where the 
turbulence field is more uniform, should show a corresponding increase in sta- 
tionarity. Such assumptions are borne out by experimental evidence. Chart re- 
cordings of image motion and image intensity over long periods of time show 
little apparent change in degree of intensity fluctuations while image motion is 
observed to periodically increase and then decrease. Measurements of finite 
record length of these quantities are nevertheless true measurements. Unfor- 
tunately, they may bear no relation to measurements of the same quantities taken 
a few minutes later. Since assignment of error bars to nonstationary data points 
would be subjective at best, let it suffice to say that the data reported here should 
be considered typical rather than average and more representative than exact. 
Noticeable spread in plots of data values, especially for image motion may be due 
largely to this nonstationarity. 
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standard statistical techniques provide a means to estimate the effect of 
finite data sample length. Data points obtained from the Stellar Image Monitor 
apparatus are random variables. A succession of these values represent a sample 
function of the stochastic processes of image intensity, motion, and size. The 
computed statistics are sample means, variances, and higher moments. For 
a number, N, of independent observations of a normally distributed random 
variable, x, it is known that the sample variance obeys a chi-square distribution 
with n degrees of freedom , i . e . , 


where 


and 



n = N - 1, 

= sample variance , 
cr^ = actual variance , 

= chi-square random variable. 


(82) 


Stellar image position at a particular time is an example of the normally distri- 
buted random variable. For weak irradiance fluctuations, as observed in the SIM 
data, image intensity can also be approximated by a normal distribution. Statis- 
tics of chi-square determine the 1 - a confidence interval for sample variance, 
S^, such that 

s2 . (83) 

2 ^2 
Xn; a/ 2 Xp ; l-a/2 

That is, the actual variance is known to fall within the above interval with a 
confidence of 100 (1 - a) percent. The key parameter that determines the above 
interval is N, the number of independent samples. In the SIM data reduction, a 
sample is taken every millisecond. Only some of these are independent, however. 
N is given by the expression 
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( 84 ) 


N = 2B T , 

eq ’ 

where 

B = equivalent white noise bandwidth 

and 

T = record length . 

In turn,iB^^ | can be estimated from the relation 


where 


B 


eq 


Rx(0) 

2 (* (r) dr 


R (r) = autocorrelation for random variable . 

X ' 


(85) 


Typical values of R^('t") for image intensity and motion lead to estimates for B^^ 
and in turn N for a given sample record length, T. The value of chi-square for a 
particular N and a is obtained from standard tables. Results are given below 
for a 90 -percent confidence interval for each signal and a sample record length 
of 20 seconds. 

Table 6 



It is apparent that from statistical considerations alone, a large range in image 
motion variance should be expected. On the other hand, the true irradiance vari- 
ance is located in a narrow range (± 5 percent) of the measured value. 
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In modulation index statistics, the important parameter is mean value. In a 
similar manner a 90-percent confidence interval can be constructed around the true 
value of modulation index mean. In this case, the statistic of interest is the stu- 
dent t distribution in terms of random variable t^, by 


where 


and 


( X - 

■■■■ 

s 

X = sample mean, 
u = actual mean , 

X 

n = degrees of freedom = N - 1 , 

N = no. of independent samples , 
s = square root of sample variance . 


( 86 ) 


The interval becomes 


X 


^ ^n,-g/2 

/N 




< X + 


v/N 


(87) 


The number of independent samples for modulation index is approximately twice 
that for image motion or about 100 for a 20 -second record. Again consulting 
statistical tables, it is possible to construct the 90-percent confidence intervals 
in terms of x and S, so that 

[x - .2 s < < X + .2s]. 

In essentially all SIM data, S was found to be less than +10 percent of y. That is, 

[ .98 X < jji^ < 1.02 x] . (86) 

Thus from statistical considerations, modulation index mean should be attainable 
with ±2 percent accuracy. 
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Tests of the electronic detection, processing, and recording subsystems 
identified the photomultiplier tube and the tape recorder-reproduce unit as pro- 
viding the limiting signal-to-noise ratios. Noise sources in the photomultiplier 
were dark current and shot noise associated with the average stellar photocurrent. 
Dark current specifications for the EMI 9558B tube used were equivalent to 1.2 x 
10” watts input. The average stellar signal (first magnitude star) for the full 
0. 152-meter aperture and optical bandwidth of 0. 1 micron was approximately 5 x 10“ 
watts indicating no trouble with dark current. A calculation including shot noise 
indicated that for the same parameters assumed above, the average photocurrent 
signal-to-noise ratio was 44 db. Thus for this case the analog tape system speci- 
fication of 43 db was the limiting signal-to-noise ratio. If limiting signal-to-noise 
was taken as 40 db, instead, the associated noise error bars would be one percent 
For smaller apertures the shot noise contribution was significant, however. 
Reduction in aperture to 0.0245 meters results in an average signal-to-noise ratio 
of 28 db. This is the limiting ratio for smaller aperture data. 

Presentation and Discussion of Results 

Stellar data reported here was obtained on various clear evenings from 
fall 1968 to spring 1969. The observation point was the Goddard Optical Research 
Facility located approximately 4.5 miles from Goddard Space Flight Center in an 
open field. The Stellar Image Monitor apparatus was placed on a 12 -inch 
high concrete slab about 20 feet from a 20- by 20 -foot trailer. This 
building supplied the necessary power and support capability. Time and reference 
signals were available from a timing system located in an adjacent building. 

The telescope and tracking mount were not enclosed by a dome or other structure. 
Standard alignment procedures were followed to assure accurate stellar tracking. 
Listed below are the nights on which data was taken and analyzed. 
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Table 7 


Nights of Stellar Observations 


Night 

(universal time) 

Date 

(local time) 

269 

24 Sep 68 

292 

7 Oct 68 

297 

12 Oct 68 

301 

16 Oct 68 

324 

18 Nov 68 

036 

4 Feb 69 

037 

5 Feb 69 

045 

13 Feb 69 

053 

21 Feb 69 

063 

3 Mar 69 

064 

4 Mar 69 

065 

5 Mar 69 

077 

17 Mar 69 

087 

27 Mar 69 


Only representative data and summary information will be presented and 
discussed here since it would be impractical to supply everything available. As 
mentioned previously, a concentration will be placed on that data most useful in 
examining the theory. 

A chart record of a typical data run is shown in Figure 18 . Experimental 
parameters that characterize the data run are hsted in Table 8. Timing marks 
at one second intervals are located along the top. Appropriate calibrations are 
indicated for the vertical scale on each recording. The modulation index track 
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has indications of percent modulation and the image motion track is calibrated in arc 
seconds. The corresponding computer statistical record is summarized in Table 8. 

In addition to the statistics listed, information on the probability density, cumulative 
density, autocorrelation, and spectral density functions is available for each data run. 

Probability density functions for irradiance and image motion statistics of 
data run 517 are plotted in Figures 19 and 20. Notice the Gaussian shaped curve for 
image motion and the long trail on the irradiance curve, indicating a log-normal 
behavior. An accurate test for normal behavior involves plotting the cumulative 
density distributions on Gaussian probability paper. This is done in Figures 21 and 
22 where the ordinate is either log- irradiance in relative units or image motion 
in arc seconds. The abscissa is constructed so that a normal random variable 
would plot as a straight line. Both quantities exhibit good fits to normal behavior 
except for edge effects where system noise and dynamic range constraints become 
important. The slope of these curves may be used to calculate log-amplitude 
variance and image motion variance. 
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Table 8 


Statistics for Data Run 118 

GENERAL PARAMETERS: 

1. DATA RUN 118 

2. STAR: ALDEBARAN NO. 168 

3. ZENITH ANGLE (ZA) = 29.7° 

4. SEC (ZA) = 1.15 

5. LOG^ [SEC (ZA)] = .141 

6. AZIMUTH = 112° 

7. 20 SECOND DATA RECORD: 301-08-46-36-56 

8. SPATIAL FREQUENCY = 0.047 CYCLES ARC SECOND 

9. APERTURE DIAMETER = 6 INCHES 

STATISTICAL SUMMARY: 

INTENSITY 

1. CIV= 0.26 

2. o-| = 0.018 

3. v = 2.84 VOLTS 

4. LOG^ (CIV) = -1.34 

5. V =0.40 VOLTS 

B 

MODULATION INDEX TRACK 

1. a /v = 0.05 

2. V = 3.9 VOLTS 

3. Vp = 0.47 VOLTS 
IMAGE MOTION TRACK 

1. V = 4.2 VOLTS 

2. 0-2 = 0.18 (ARC SECONDS)^ 

3. cr = 0.43 ARC SECONDS 


96 




NIGHT 63 

J DATA RUN 517 

•• STAR 267 

• 5 ZENITH ANGLE 56.1® 


CIV=-pL= .43 

• o 

• • 

• • 


0 ■ 1 

0 1.0 

2.0 3.0 

IRRADIANCE 

(RELATIVE UNITS) 

Figure 19. Irradiance probability density function 
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Figure 20. image motion probability density function 
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Figure 21. IRRADIANCE cumulative density function 
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Figure 22. Image motion cumulative density function 


A summary of irradiance statistics for all observations on nights 37 and 65 
is presented in Figure 23. All data was taken with the full 0. 152-meter aperture 
and is plotted in the form log^ (CIV) versus log^ [secant (zenith angle)]. The 
straight line represents the linear least squares fit to all data points with zenith 
angles less than 60 degrees. The calculated slopes agree well with the experi- 
mental value of sec^-2 (ZA) reported by Protheroe. For larger zenith angles 
the expected saturation and data scatter effects are evident. Vertical axis inter- 
cepts yield the values of CIV for zenith viewing, 20 percent and 14 percent 
respectively on the two nights. Protheroe reported average values of 14 per- 
cent and 16 percent for summer and winter viewing with the same size aperture. 
In Figure 24 zenith -angle dependence of log-amplitude variance, cr|, is plotted 
as well for night 65. Zenith angle dependence reflects the proportionality of 
to CIV^. A typical error bar is shown for one point in each graph. The 
remaining data scatter is likely due to nonstationarity since observations en- 
compassed a period of three or four hours on each night. 
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Figure 24. Log-amplitude statistics 


Aperture dependence of irradiance statistics is reflected in the three series 
of observations in Figure 25. A series of five data runs with aperture diameters 
from 0.0508 meters to 0.152 meters was taken on each of the three successive 
nights in March. The solid line through each data set represents the best expo- 
nential least squares fit of the form 


CIV^ = A exp (- B diameter). 


The parameters of the fit are A, B, and R, a correlation coefficient representing 
the accuracy of data fit to the exponential. The parameter B has units of meters”^ 
and can be related to the irradiance correlation distance, /Xz, by the technique 
described in Chapter III. The calculated values must be corrected for zenith 
angle of observation in order to arrive at the desired quantity V ^ , 


where 

and 


/Xz “ V s e c 6^ , 

6 = zenith angle 

Hq = altitude most significant in determining irradiance fluctuations. 
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Figure 25. Irradiance aperture dependence 
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Results are tabulated below for the three series of runs. Height at which maximum 
wind velocity occurred is given, as is the height for the temperature inversion cor- 
responding to the tropopause. The meteorological parameters were obtained from bal- 
loon rawinsonde probes launched from near Dulles International Airport within a few 
hours of each series of runs. Although the launching site is about 80 kilometers 
from the point of stellar observations, the upper altitude parameters can be con- 
sidered fairly typical for the entire Washington, D.C. , vicinity. 


There is no clear correlation between measured height, Hq, and heights at 
which metero logical evidence for turbulence exists. The statement about correla- 
tion distance was that it should be of the order of » not necessarily equal to it, 
however. All heights are within a factor of two of each other, indicating that effects 
located at these higher altitudes are contributing most significantly to irradiance 
statistics. 
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Table 9 


Irradiance Correlation Distance and Meterological Parameters 



Exponential Fit 

Meterological Parameters 

Night 

B 

/Xz 

Am, 

Ho 

height of 
tropopause 
(from temperature) 
data 

height of 
max. wind 
velocity 

max. 

wind 

velocity 


(meters"* ) 

(meters) 

(meters) 

(kilometers) 

(kilometers) 

(kilometers) 

(meters/second) 

63 

9.91 

0.134 

0.096 

18.4 

10.1 

9,2 

40 

64 

11.9 

0,112 

0.083 

13.8 

10.0 

14.5 

18 

65 

8.54 

0,156 

0.106 

22.5 

11.6 

10.8 

66 





Frequency domain behavior of irradiance is reflected in Figures 26 and 27 
and also in the figures of Chapter V. The latter are of particular interest since 
they show spectrum shape for near -zenith viewing on several different nights. 

They will be discussed in connection with a model for strength of turbulence in 
the latter chapter. It has been shown that frequency spectrum of log-amplitude 
is nearly identical to that for irradiance under the conditions of small For 
the starlight observations reported here this is certainly the case. Essentially 
all spectral densities of irradiance have the same type of shape. They tend to 
increase slowly from zero frequency to a broad peak at about 10 Hz to 20 Hz. 
Beyond that they obey a near exponential decrease to the system noise level 
(40 db) at a few hundred Hz. It is reasonable to assume that turbulence elements 
of size, /Vz, which are most responsible for the irradiance fluctuations also 
are responsible for the peak. Without prior knowledge of the proportionality 
constant relating frequency to />Tz, the height z cannot be determined here. As 
zenith angle increases, the peak moves toward the origin and bandwidth of the 
entire spectrum is reduced, roughly in proportion to 1/ sec ^ . As aperture 
diameter decreases the peak moves to higher frequencies and bandwidth increases. 
Samples of zenith angle and aperture dependence are shown in Figures 26 and 27 
respectively. These spectral densities are normalized to unit variance for ease 
of comparison. 

An exponential of the form 

CIV^ = A exp (- B f req. ) 

can be fitted to the data beyond the low frequency peak and B ^ used as a measure 
of bandwidth. When this is done the ratio of bandwidths for the different zenith 
angles of Figure 26 is found to be 1.8. This compares with 2.0 for the ratio of 
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Figure 26. Example of irradiance spectral density zenith angle dependence 



Figure 27. Example of irradiance spectral density aperture dependence 
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square root secant zenith angle for the two stars. Aperture reduction in band- 
width is about 2:1 in passing from the 0.0245-meter aperture diameter to the full 
0.152 meter. No strict power law dependence is expected and this matter will be 
taken up again in the next chapter. Relationship of bandwidth, B ^ , and upper 
altitude wind velocity is dramatically demonstrated in the following table. Listed 
here are exponential fit parameter, B ^ , and R, correlation coefficient, along with 
the maximum wind velocity for various nights. All spectral densities were com- 
puted for near-zenith observations with the 0.152-meter aperture. This prevented 
confusion between zenith angle and aperture dependence. 


Table 10 

Irradiance Bandwidth and Upper Altitude Wind Velocity 


1 

Night 

Data 

Run 

Zenith 

Angle 

CIV 

Exponential Fit Parameters 

\ (max) 
(meters/second) 

A 

B 

R 

1/-B 

292 

61 

6 

.085 

.0135 

-.0229 

-.930 

43.7 

21 

297 

87 

7 

.112 

.0204 

-.0187 

-.982 

53.5 

38 

301 

121 

14 

.244 

.0338 

-.0273 

-.978 

36.6 

32 

036 

424 

10 

.252 

.0102 

-.0186 

-.931 

53.8 

54 

037 

449 

7 

.220 

.0354 

-.0263 

-.986 

38.0 

51 

063 

524 

18 

.199 

.0205 

-.0171 

-.985 

58.5 

40 

064 

547 

11 

.139 

.0317 

-.0323 

-.991 

30.96 

18 

065 

584 

18 

.135 

.0316 

-.0149 

-.957 

67.11 

66 

077 

623 

14 

.242 

.0921 

-.0118 

-.962 

84.5 

97 

087 

633 

25 

.180 

.0257 

-.0215 

-.978 

46.5 

38 


While irradiance statistics were found to be fairly similar from night to 
night, image motion and image size statistics changed markedly within even a 
few hours on a particular night. This was to be expected due to the higher degree 




of nonstationarity of the lower atmosphere, where these quantities are primarily 
determined. Increase in statistical errors for the measurement of motion and 
size also played a large role. Two sets of image motion statistics for night 37 
are shown in Figure 28. The first set was taken just after dusk and the second 
set an hour later. Note the factor-of-two reduction in for all zenith angles 
during the later period. Lines drawn through the data represent the best fitting 
line of the form 


0-2 = M secant (zenith angle). 

Value of the slope m was 1.4 and 0.65 (arc seconds)^ for the first and second sets 
of observations respectively. M is the value of for zenith viewing. For 
essentially all observations on all nights the range of M was 0,4 to 1.6 (arc seconds)^. 
Statistical error bars drawn for a typical point in Figure 28 indicate that nearly 
all data scatter can be accounted for on this basis alone. Before plotting, each 
data value of cr^ was corrected for a noise component most probably due to 
chopper pattern centering error. Fortunately the noise was periodic (8 Hz) and 
its magnitude could be accurately measured from the spectral density of the 
voltage output proportional to image motion. Magnitude of the noise was generally 
such that it contributed 5 to 20 percent of image motion variance. All values of 
(7^ and spectral densities plotted here have been corrected individually for the 
noise component. 

Data on aperture dependence of image motion was taken simultaneously with 
the three sets of data on irradiance aperture dependence discussed previously. 

Each data set was taken over a 15 -minute period and consisted of five different 
aperture diameters evenly spaced from 0.051 to 0.152 meters. As shown in Figure 
29 a least squares power curve of the form 
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IMAGE MOTION VARIANCE CT^ (ARC SECONDS) 


1.4 SEC (ZA) 


^90% CONFIDENCE INTERVAL 


0.65 SEC (ZA) 


TIME PERIOD 

OHRS 19 MIN TO 

1 HR 21 MIN 

2 HRS 13 MIN TO 
3 HRS 9 MIN 


ZENITH ANGLE 


20*^ 30° . 40° 


SECANT (ZENITH ANGLE) 
Figure 28- Image motion statistics 


LEAST SQUARES FIT 
o NIGHT 64 crm = 2.7D°'’®R=-0,77 

• NIGHT 63 CTm= 1.2 ='0.30 

a NIGHT 65 CTm =0.86D°^^R=-0.68 


90% CONFIDENCE INTERVAL 



APERTURE DIAMETER D (METERS) 
Figure 29. Image motion aperture dependence 






O' 2 = A (diameter)® 

was fit to each data set. The calculated power dependence, B, compares roughly 
in each set with the theoretical value of -1/3. Correlation coefficient, R, is 
low, especially for night 63 data, indicating data scatter. Once again, most 
scatter can be explained by statistical error bars based on the length of each sample 
and the number of independent samples. 

The frequency spectrum of image motion for data taken with a 0. 152-meter 
diameter aperture is characterized by a rapid exponential decrease from zero 
frequency. Essentially all motion is confined to frequencies below 20 Hz. A 
typical example is shown in Figure 30 and several more are included in the 
Figures of Chapter V. A decrease in aperture or zenith angle tends to flatten 
the spectrum and extend it to somewhat higher frequencies. This data and its 
relation to strength of turbulence analysis will be examined in greater detail in 
the next chapter. 
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Image size information is derived from the measured modulus (MTF) of an 
optical transfer function. Plotted in Figure 31 are measured MTF values for two 
sequences of stellar observations on night 65. Each plotted value is determined 
by the calibration of modulation-index voltage. Observations were taken on each 
star at five spatial frequencies during a total time period of 15 minutes. 

A smooth curve was then fitted through the data points and the modulus = 1.0 point 
at zero spatial frequency. The resultant MTF does not yield information on image 
fine structure but provides an estimate of the average short-term profile. The 
profile plotted in Figure 32 is the computer-calculated Fourier transform of 
the observed MTF for star 193 on night 65, corrected for optical system MTF 
by division. The calculated profile for star 257 on the same night was not plotted 
since it fit almost exactly the star-193 result. Observed modulus values for the 
two stars show a distinct difference (10 to 20 percent) for intermediate spatial 
frequencies, however. This indicates that a 10 percent error in modulus value would 
not be reflected in a correspondingly large error in image size. The Airy dif- 
fraction pattern for an f/8, 0.1 52 -me ter diameter, circular objective in 0.5 -micron 
light is plotted in Figure 32 for comparison. Width of both profiles is measured 
at the l/e points. The Airy pattern was generated from the diffraction limited 
MTF of Figure 17 by the same computer routine used on stellar data. This serves 
to test the computational approximations used. Note the expected ring structure 
on the Airy pattern. Ring structure and negative intensities on the stellar profile 
are caused by truncation of the Fourier Transform at 0.6 cycles per arc second. 

This is where data on modulus ended. This truncation does not sensibly effect 
image width, however. 


I 
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An estimate of long-term profile is obtained as discussed in the previous 
chapter. For the long term, the key parameter is image motion variance, cr^ . 
Two profiles corresponding to the data sets of star 193 and star 257 from night 65 
are plotted in Figure 33 . The large difference in profile width reflects the differ- 
ence in cr^ values for the two stars. Note the disappearance of ring-like struc- 
ture. This disappearance can be attributed to the smoothing action of the transfer 
function computed from image motion variance. Variance of image profile, cr^ , 
can here be safely used as a measure of width since there are no false contribu- 
tions from ring-like structure. Departure of ratio of the two cr^ 's from 
the ratio of secant (zenith angle) is expected since ratio shows this departure. 
If each profile were strictly Gaussian, its width at the l/e points would be 
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ANGULAR DEVIATION (ARC SECONDS) 
Figure 33. Long-term image profiles 
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Table 11 


Stellar Image Sizes 



This behavior is also roughly confirmed. A summary of all image-size results 
is given in Table 11. These results reflect the wide variety of conditions present 

on the various observing nights. Note again the large dependence of the long- 
term image size on strength of image motion, since the method of obtaining 
or^ involves essentially addition of cr^ and the short term o-^. 
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CHAPTER V 


DETERMINATION OF TURBULENT PROFILE 
Best Present Estimates of Profile 

The common factor for all of the theoretical analysis reported in Chapter III 
is an integral expression containing C^(h), the refractive index structure constant, 
as a function of altitude. Lack of both experimental and theoretical knowledge 
of this profile is the single greatest obstacle to understanding vertical propagation. 
If the profile were known, optical propagation theory could be tested in any number 
of interesting experiments. Yet what data there is suggests that no one series of 
measurements to date has been sufficient to accurately describe the profile. For, 
as soon as the measurements were taken, the atmosphere changed. This is the 
old problem of nonstationarity of the turbulent process. Lawrence, Ochs, and 
Clifford"^^ in a series of nighttime measurements with an airborne temperature 
sensor noted several orders of magnitude fluctuation in C^(h) for essentially the 
same altitude. All these measurements were taken within four kilometers of ground 
level and probably reflect the presence of time-var3dng local disturbances. Turbu- 
lence at higher levels should be relatively more stationary. Of course, the average 
behavior of atmospheric turbulence could be described by a careful series of 
thermal and optical measurements as a function of altitude. A program to do this 
would be difficult and costly and would not always generate results when needed. 
Alternately, if the theoretical dependence on C^(h) is accepted as a working hypo- 
thesis, much can be learned from ground-based observations. A useful concept in 
conjunction with these observations is a turbulent profile model, i.e. a functional 
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dependence for C^(h). Model parameters can then be adjusted for correct predic- 
tion of ground-based stellar data. In particular, simultaneous Stellar Image Monitor 
observations of both irradiance and phase statistics should permit a more accurate 
determination of model parameters than previous data. Neither direct measure- 
ment nor ground-based observations solve the real problem of nonstationarity, 
but the technique of modeling should make it possible to describe average or at 
least typical behavior of the turbulent profile. 

Any proposed model should reflect probable meteorological factors at the 
origin of turbulence. As noted previously, essentially all refractive -index varia- 
tions are accounted for by small-scale temperature fluctuations. These in turn 
can be related to the pertinent larger scale meteorological parameters. The intent 
is not to construct an actual profile from meteorological parameters, but merely 
to show their influence. Equation (13) of Chapter II related C^(h) to G^(h). This 
can be approximated for visible wavelengths by 


C2 (h) ^ 10-12 (^\ Cl (h). 


(89) 


In turn, C^(h) can be expressed from turbulence theory 
meteorological parameters for a particular altitude by 

C2 Ch) 3 A2 


19,32 


in terms of gross 


^2 'S y2 


(90) 


where 


e = rate of energy per unit mass dissipated by viscous friction, 
y = average vertical gradient of potential temperature, 6 


-<U>' 

P = average vertical wind shear, 

and 

= empirical constant. 

Numerous measurements of e are available for at least the lower altitudes. 

2 3 

Typical values, as quoted by Hufnagel^'*, show a decrease from 300 cm /second 
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at ground level to 0.07 cm^ /second^ at 10 kilometers and then an abrupt rise, 
although data is scarce and quite variable in this region. Potential temperature, 0 , 
is defined by 

6 ^ T ~ rh, 

in terms of actual temperature andP , the adiabatic lapse rate. For dry air 


r = ~ 9.8 ^ C /kilome t e r . 


The relationship between 3T/3h and T is a measure of stability, such that 

^ T ^ P ^ 1 

7 — < 1 unstable , 
o h 


B T 
B h 


= 0 


neutral, 


and 


> r s table, 
B h 


— - > 0 inversion, 

B h 


Near neutral and unstable conditions favor the development of turbulence. Values 

44 

of /3 are the least well known. As Hulett pointed out, small-scale shears between 
layers a few tens of meters apart tend not to be included in the reported data. 

Shears of this dimension are probably the most effective in controlling optical 
effects. The three parameters of equation (90) are height dependent and also a 
function of one another. Thus it is difficult to isolate a single cause of turbulence. 
The important information is that the most probable mechanism for turbulence 
generation, at least in the lower altitudes, is vertical temperature gradient, vertical 
wind shear, or some combination of the two. 


Upper -altitude investigation (above several kilometers) of aircraft size turbu- 
lence by Clodman, Morgan, and Ball"^^ suggests that wind shear can also combine 
with a gravity -wave mechanism to yield enhanced turbulence. Gravity waves re- 
sult from large air mass flow over obstacles such as mountain ranges or even 
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small hills. They are only imperfectly understood at present. It is known that 
as these waves break down and interact with each other and the main flow of wind 
itself, severe turbulence may result. These waves tend to propagate to regions 
of lower density, meaning higher altitudes. During this propagation, gravity wave 
amplitude tends to build up and at the tropopause quite a strong disturbance 
could result from interaction with shear -generated turbulence. The tropopause, 
the region of temperature inversion at approximately the 200 -millibar pressure 
level, has historically been considered a source of turbulence. Whether it acts 
as a source or permits amplification of existing structure is still unknown. 

Perhaps the most significant feature of this region is the existance of jet stream 
flows. Velocities of maximum winds in this flow range from 20 to over 100 meters 
per second. The jet stream vicinity is a region of strong horizontal and vertical 
wind shears. This may present the source of turbulence energy input which cascades 
down to the smaller scale phenomena affecting optical propagation. Observations'^^ 
indicate that the maximum turbulence occurs approximately ±40 millibars from 
the jet stream core. This roughly corresponds to 1 kilometer below and 2 kilom- 
eters above the core. Conditions for turbulence generation at higher altitudes are 
not very sensitive to diurnal variations, as are conditions in the surface boundary 
layer. The temperature field is much more stable near the ground at night than 
during the day and strength of turbulence is often reduced by a factor of 100. 

A typical example of vertical profiles of temperature and wind velocity are 
shown in Figure 46. This is data from a rawinsonde observation for March 6, 1969. 
The point of observation was, Sterling, Virginia, and time release was 6:15 pm EST. 
The midpoint of the flight corresponded to the beginning of stellar observations on 
that particular night. No attempt at smoothing the reported profile was made 
The temperature profile indicates a small inversion layer at the surface and then 
an abrupt switch to near by neutral conditions. Stability with little change in the 
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value of temperature gradient is evident in the tropopause at 11.6 kilometers, 
except for a disturbed temperature layer near 1 kilometer. Above the tropopause 
the gross temperature profile exhibits fluctuations in the value of temperature 
gradient. Wind profile maximizes at 10,6 kilometers and a speed of 66 meters 
per second. This is well-developed jet stream flow and represents almost a linear 
increase from the surface winds of two meters per second. Above the tropopause 
wind velocity oscillates in magnitude indicating the presence of strong vertical 
wind shear. Wind direction remains fairly constant for all levels. Data points 
above the tropopause are separated by a kilometer or more and thus reveal 
little about fine structure. Smaller scale wind shears and temperature gradients 
are most probably associated with these regions of large-scale shears and gradients, 
however. 

There are several important points to be derived from the meteorological 
theory and data. 

(1) Little meteorological evidence exists for a thin, highly turbulent layer 
responsible for all the observed optical effects. 

(2) Mechanisms do exist at nearly all levels of the atmosphere for the de- 
velopment of turbulence. The combined factors of wind shear, potential 
temperature gradient, and gravity waves are the most likely sources. 

(3) A turbulence boundary layer exists near the earth^s surface in the first 
tens to hundreds of meters. 

(4) Intermediate altitudes (below the tropopause) are thermally stable in 
general (especially at night) and less turbulence is expected. 

(5) The region of the tropopause and beyond is characterized by conditions 
of possible instability and a high velocity jet stream flow. Turbulence 
is expected to increase in this region, yet not surpass the maximum 
near -ground values. 
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Solutions for Model Parameters 

It is desirable to have a functional dependence on altitude for strength of 
turbulence, C^(h), in order to facilitate comparison of data with theory* At the 
outset it should be realized that no model could possibly account for the fine 
structure and all the disturbed layers that are thought to exist. The model should, 
however, take into account the gross meteorological factors mentioned above and 
conform to the average or at least typical behavior of the atmosphere. As a first 
approach, consider a model that uses an exponential decrease with altitude to 
describe the envelope of C^(h) fluctuations. Re tractive -index varies linearly with 
density which in turn is proportional to exp (- h/7) for a standard atmosphere. 
C^(h) goes as density squared, thus could be modeled by 


where 

and 


C2 (h) exp (-h/3.5) , 

= ground strength of turbulence 
h = altitude in kilometers . 


(91) 


The hypothesis is that turbulence due to essentially the same mechanism as that 
near the ground would be damped out from lack of material at higher altitudes. 

This certainly does not represent the full behavior of C^(h) and appears "^"^to over- 
estimate the contribution from intermediate levels below the tropopause. A more 
accurate model in the sense of agreement with data can be constructed from the 
curve supplied by Hufnagel, as shown in Figure 34. His values for a clear evening 
arise from meteorological relations and stellar observations. The curve is charac- 
terized by an exponential type decrease in the lower atmosphere, a gradual leveling 
off, and a sharp spike at the tropopause. A mathematical approximation to this 
curve can be written 


exp C - h/H ) 

C2 (h) = A ^ + B S (h - , 


(92) 
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Illlllllllll 


where 

A = 9.1 X 10"^^ meters , 

H =1 kilometer , 

0 

B - 2 X 10“^^meters”^^^ , 

, H, = 12 kilometers , 

and 1 

a = 5/3 . 

The model of high turbulence at the tropopause represented here by the delta 
function, was included by Hufnagel to explain stellar irradiance fluctuations. The 
value of and the power of altitude in the denominator indicate a more rapid 
roll-off than the atmospheric density exponential. The density exponential is 
plotted in Figure 34 for comparison. Hufnagel also presents daytime data which 
more closely corresponds to the exponential density decrease. There are some 
scattered data points in the Literature that could be added to Figure 34 for com- 
parison. However, these results were most often obtained by inaccurate means 
and invariably were daytime results. To the author ^s knowledge one set of useful night- 
time results does exsist, at least for the lower atmosphere. This is the set of high- 
speed temperature measurements made from an aircraft and reported by Lawrence, 
Ochs, and Clifford"^^ . These results seem fairly dependable since system noise 
level was less than 10'^® meters"^ ^ . The average of their measurements agrees 
well with the Hufnagel model although individual C^(h) values fluctuate orders of 
magnitude around this average. The net result of the Hufnagel model is a division 
of C^(h) profile into low-altitude and high-altitude structure and an estimate of the 
relative contribution of each altitude regime. Hufnagel’ s quoted result for the 
integral over the whole profile is 4.3 x in units of meters . 

The relationship between the two regions of the Hufnagel model and a check 
on self consistency of the theory can be made by inserting the model into equations 
developed in Chapter HI. All the theoretical expressions for optical statistics 
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depend on an integral over C^(h) , the refractive -index structure constant. The 
form of the integral varies with the statistic being considered. In the case of 
log-amplitude or irradiance, the integral 




1 “ 




hS /6 ^ 


is present and indicates a non-uniform weighting of C^(h) values. The largest 
contributions to the integral are at higher altitudes where the product C^(h) h^^® 
takes on its largest values. Image motion, modulation transfer function, and image 
size statistics contain the integral 


^ CD 

C2 (h) dh, 

which, of course, takes on its largest values where C2(h) is the largest, presumably 
near the ground. 

Assuming a model of the form fitted to the Hufnagel profile, an estimate of the 
parameters of this model can be obtained from the theoretical expressions for 
log-amplitude and image motion. Of these parameters the ground strength of 
turbulence is the best substantiated at approximately 5 x meters . 

Even this value can exhibit a large variation with locality and time of measure- 
ment. It is generally agreed from various optical and thermal measurements 
that this value is accurate ± one order of magnitude near the ground at night. The 
height of the tropopause can also be determined with some assurance from 
meteorological measurements on a given night. The range is approximately 9 to 
14 kilometers. The two remaining parameters, and B, can be arrived at through 
the following process. First the anal 3 d;ical expression is inserted in the equa- 
tions for log-amplitude variance and image motion variance and the integrals are 
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solved. Then the measured values of and cr^ and the assumed values for A 
and are inserted. The result is two equations in the two unknowns, Hq and B, 
which can be readily solved. Thus a model of the turbulent profile can be con- 
structed for a particular time on a particular night. This model can be checked 
for self consistency by insertion into the equation for image size. For ease of 
mathematical treatment the power of the denominator in equation (92) is changed to 
2/3. This affects the fit of the model to HufnageFs curve but still gives order of 
magnitude results. The integrals are now readily evaluated in terms of gamma 
functions, as follows: 


m 1 


C2 (h) d h , 


where 


and 


where 


Then, 


Kj = .56D‘1'3 SEC(ZA), 


^1 = ^2 ( h ) 


.09 SEC2-4 (Z A). 


exp(-h^H.) r 

= K, A 1 . d h + X, BS(h-H,)dh 

J h2/3 J 

0 

exp (- z) d z + J B 


(93) 


(94) 


= [ AH^i '3 r(l'3) +B], 


and 


exp(-h/HJ f“ 

? - - A hS/6 d h + /<2 B S (h - Hj) hS d h 

Jn J 


- '"2 
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exp (- z) z^/^ d z + /<2 


= K, 


AH^/6 


f 

• 4 ) 


= /<2 [AH7/6 r (7/6) + B HS/6] , 

where 


I 


03 


exp (- z) z^~ ^ d z - r (4>) • 


(95) 


As a typical example, consider data run 584 from night 65, for which 


and 


diameter, D= .152 meters, 
SEC (Z A) z: 1.06, 

X. = .5 X 10”® meters, 

- 11.6 kilometers, 

A = 4 X 10”^^ meters”^'^^* 


Noting that F (7/6) = 0.923 and F (1/3) = 2.68, the two equations can be solved to 
obtain 


and 


H^ = 112 meters 
B = 6.7x 10”^® meters“^^^. 


Equation (76) for image size under the assumption of a Gaussian beam profile 
is also expressible in terms of these parameters, as follows: 


where 


K 


3 


a 


2 

S 



(h) d h , 


. 98 SEC (Z A) , 


(96) 
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and 


o'i = 


A, 

n 


exp (- h/Hp) 
h27i 


I 


^ k\ B S (h -Hj) d h 


= /<3 [AHJ/^ r(l/3) +B]. 

Inserting the above calculated values for Hq and B yields 

0-2 = 2-9 (arc seconds)^. 

This compares with the value of o-^ = 1,1 (arc seconds)^ computed from long- 
term MTF obtained nearly simultaneously with a-2 ando*^. The values of 
and B also compare roughly with the Hufnagel fit values of one kilometer and 
2 X 10”^^ meters”^' ^ . Estimates of Cisj^(h) profile obtained in this fashion are, of 
course, restricted by the applicability of the model. 

Spectral densities of irradiance and image motion provide a more sensitive 
technique to arrive at an estimate of specific shape of C^(h) profile. As noted in 
Chapter III both spectral densities are expressible in terms of an integral of the 
form 


c2 (h) g (h), h, D, e, f} d h, (9'^) 

where 3 { • } is a nonlinear function dependent on some or all of the following 
parameters: 

(1) Altitude, h. 

(2) Wind velocity profile, (h) . 

(3) Telescope aperture diameter, D. 

(4) Stellar zenith angle , 6- 

(5) Temporal frequency under consideration, f. 
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The 9 {• } function in both irradiance and image motion statistics is a fairly 
complicated integral expression. Peskoff and Fried have solved a simpler 
version (no frequency dependence) of the irradiance equation. The result in the 
form of a hypergeometric function is useful only when expressed in a power series 
approximation and is restricted to the point detector case. Another approach as 
used in this thesis is numerical integration. The 3 { • } equation in each case for a 
finite aperture and explicit dependence on wind velocity becomes a tabulated function. 

Reasonably accurate values of wind velocity are available from the Weather Bureau 
and numerical integration is done by digital computer. The result is a linear 
integral equation connecting the measured values of or ClV^to the unknown 
profile, (h). If the integrations are replaced by a summation over a number of 
finite h values, the result is a set of linear simultaneous equations in the unknown 
(h) values. In principle this set could be solved for the desired unknowns by 
standard mathematical techniques. This approach was tried and very anomalous 
oscillating solutions were obtained. Apparently there is a problem with data signal- 
to-noise ratio and, to a lesser extent, with uniqueness of solution. In any event, 
considering the error bars on spectral density estimates and the wind velocity 
profiles, it seemed more appropriate to limit the number of degrees of freedom 
and achieve a rougher estimate of the profile. This was done by using a model for 
(h) based on the Hufnagel curve but containing six unknown parameters 


exp (- h 

(h) = A + B exp 

h“ 


(h 

2 j82 


( 98 ) 


The first term is the familiar combination of exponential roll-off and power of 
altitude in the denominator. This controls the low-altitude dependence. The second 
term is inserted to account for the presence of upper altitude disturbances . It is 
a Gaussian shape with magnitude, B, rms size, jB , and centered at height, Hj . The 
Gaussian shape was chosen arbitrarily as a convenient way of varying the magnitude. 
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size, and position of the representation of upper altitude effects. Another curve 
or distribution might work just as well and may even be a truer representation 
of the physical process. In the absence of prior knowledge of the actual physical 
process, one can well use the Gaussian curve. 

Presentation of Data and Construction of Several Profiles 

Realizing that image motion effects are due primarily to the strong turbulence 
of the earth^s boundary layer, the parameters, A, H^, and a were determined 
through image motion spectral density data. The three parameters were allowed 
to take on discrete values and a test spectral density was constructed by numerical 
integration of equation (62). Numerical integration consisted of summing contribu- 
tions to spectral density from discrete altitude intervals. Each test spectral 
density was evaluated at a number of discrete frequency values and subtracted 
from the corresponding data value at that frequency. Absolute value of this dif- 
ference was used as the measure of accuracy of fit to the data. For image motion, 
altitude intervals of 10 meters were used in the first 100 meters and then intervals 
of 100 meters were used to an altitude of one kilometer. The interval range was 
then broadened to 500 meters and the maximum altitude examined was 9. 5 kilometers. 

Contributions to spectral density were examined every 0.5 Hz from zero frequency 
to 15 Hz. A test spectral density was computed for the C^(h) profile resulting 
for every possible combination of the three parameters. In a typical case each 
parameter was allowed to take on 10 different values. This meant that the results of 
1000 combinations were tested. Both the numerical integrations and testing were 
conducted on a digital computer. The computer was programmed to select the best 
fit to spectral density data, print it out and also print out the C^(h) profile corre- 
sponding to this best fit. Also printed out were the sum of absolute differences for 
all frequencies for each test case. This enabled an assessment of fit sensitivity 
to parameter changes. Irradiance spectral density data was used in the same 
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fitting routine to determine the Gaussian shape parameters B, and ^ . Since 
irradiance data was expected to be somewhat more influenced by lower altitude 
effects than was image motion data by higher altitude effects , the image motion 
fitting was done first and the best fit parameters used in equation (98) for the full 
C^(h) profile. Then the irradiance fit was performed to find the remaining three 
parameters. Altitude intervals for the irradiance data were taken at 100 meters 
for the first kilometer and then increased to one kilometer up to a maximum 
altitude of 25 kilometers. This represented the sensible limits of the atmosphere 
although the formal integration extends from zero to infinity. Contributions to 
spectral density were examined every 10 Hz from 10 Hz to 250 Hz. At the upper 
end of the frequency scale irradiance spectral density was usually down at least 
a factor of 100 from the 10 Hz value. 

At first each parameter value was allowed to vary over a wide range. The 
ranges were than narrowed as it became obvious which values gave the best fits. 
It was found that a pure exponential decrease, while giving a close fit to image 
motion data, could not explain the shape of irradiance data. The image motion 
data required an exponential roll-off typically less than one kilometer while 
irradiance data favored roll-offs of tens of kilometers without yielding a good 
shape for even these high values. The implication was clear. Significant contri- 
butions to C^(h) were needed at higher altitudes to correctly predict irradiance 
spectral density. With further analysis it was found that the Gaussian "bump^' 
fit the data quite well. Data and the corresponding fits for six different nights 
are presented and discussed below. These nights represent a variety of wind 
conditions and presumably different distributions of turbulence. 
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The use of the 'sum of absolute differences between test and data spectral 
densities was justified on the basis of the desire to obtain a best fit to the shape 
of the curve. Absolute difference was sensitive to the lower frequencies where 
spectral density was the largest. This was somewhat desirable since here the 
signal-to~noise ratio of data was greatest. Statistical error bars were determined 
in Chapter IV. For irradiance these were ±5 percent and for image motion, ±50 
percent. These of course apply to total variance, the integral over the entire 
spectral density. A relative idea of error distribution with frequency can be ob- 
tained by comparison of two consecutive data records, as in Figures 39 and 40. 

The irradiance spectra exhibit differences of about 10 percent near 10 Hz. Although 
at frequencies near 200 Hz the relative difference is as high as 100 percent in 
places, the general shapes of the two curves match well. The two records of 
image motion show large differences (a factor of four) at the origin, but beyond 
one Hz, have nearly the same values as well as shape. Discussion of this low- 
frequency variability in terms of atmospheric, instrumental, and resolution 
effects is included in Appendix B. A mean square difference statistic was also con- 
sidered but rejected for its heavy weighting of the lower frequencies. Since a lot 
of shape information for both types of data is contained in the higher frequency 
roll-off, it was felt that absolute difference was a reasonable compromise between 
equal weighting at all frequencies and the mean square approach. 

Wind velocity information obtained from weather bureau rawinsonde flights 
was supplied in terms of speed and direction for a given altitude above the surface. 
The reported speed was that component parallel to the earth^s surface. For 
spectral density analysis using Taylor's hypothesis, the required speed component 
is that perpendicular to the optical line of sight, . For horizontal winds of 
speed Vq there will obviously be a reduction in for stellar observations at 
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appreciable zenith angles. As Young has pointed out there is also a possible 
increase in dependent on relationship of wind direction to azimuth of the stellar 
observation, by the expression 

[1 + TAN2 (Z A) SIN2 (A Z - 1/2. (99) 

where 

= wind component perpendicular to line of sight, 

Vq = horizontal wind speed , 

Z A = stellar zenith angle, 

A Z = stellar azimuth angle, 

and 

6 ^ = direction from which wind is blowing. 

Division by secant (ZA) would be the correction factor if only zenith angle effects 
were considered. The term in brackets represents the interaction of zenith angle, 
azimuth angle, and natural wind direction. This term is maximum when (AZ - 0 ^ ) 

= 90°. For this case and a zenith angle of 60° the azimuth correction essentially 
cancels the zenith angle correction. To eliminate these effects from the spectral 
density analysis only observations taken within a few degrees of zenith were used. 
For this situation, within a few percent regardless of wind direction. 

Speed, Vq , was itself determined to this same order of accuracy. Although the 
point of rawinsonde release was some 80 kilometers from the point of stellar 
observations, the release time and prevailing winds were such that the balloon 
was essentially overhead for most of the stellar observations. In addition, the 
upper altitude winds showed little change over a few hours or within a few hundred 
kilometers of the measurement site. Rawinsonde data for March 6, 1969, taken at 
stations as diverse as Wallops Island, Virginia; Richmond, Virginia; and Phila- 
delphia, Pennsylvania, and for time separations of up to 6 hours showed correla- 
tion within 20 percent with V^ values obtained at the regular Sterling, Virginia, 
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station. Below 1 kilometer altitude, correlation began to worsen indicating the 
presence of local terrain effects. Below 150 meters rawinsonde wind velocity 
was not available. Unfortunately this is the region most effective in determining 
image motion spectral density, A reasonable approximation to actual wind pro- 
file was obtained by linear interpolation from a recorded or estimated surface 
wind value to the first rawinsonde value at 150 meters. All surface wind values 
at the point of observation were less than three meters per second, for beyond 
this value telescope vibration became a serious contribution to image motion and 
data was not taken. 

The six nights for which turbulent profile models were constructed are listed 
in Table 12 along with the particular data run parameters used on each night. 
Image motion and irradiance spectral densities from these data runs are plotted 
in Figures 35 through 41. Ordinate values are o-^/Kz and CIV Vhz respectively. 
Thus the strength of turbulence induced fluctuations at a particular frequency can 
be compared from one night to the next. Corresponding wind velocity profiles 
for each night appear in Figures 42 through 47. Note the significant range of 
maximum wind speeds represented here. The continuous solid line on the spectral 
density plots is a smooth curve representation of the best computer generated fit 
to the data points. Turbulence profiles which resulted in these best fits are shown 
in Figures 48 through 53. 

Computer generated fits to image motion spectral density followed the average 
behavior of the data. In particular the sharp increase below one Hz was well 
approximated. A disturbing feature of the fit, however, was the one to two Hz 
oscillation observed in all curves. Oscillation resulted from a truncation of the 
Fourier transform from the computed autocorrelation function. The maximum lag 
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Figure 36. Spectral densities night 301 
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Figure 37, Spectral densities night 36 
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Table 12 


Data Runs for Calculation of Turbulent Profiles 


Night 

Data 

Run 

Number 

Stellar 

Zenith 

Angle 

Time 

CIV^ 

C 7 ^ 

m 

(hrs) 

(mins) 

297 

87 

70 

23 

56 

.0127 

.13 

301 

! 

121 

0 

1 — 1 

09 

05 

.0596 

1 

.37 

36 

424 

0 

0 

00 

45 

,0635 

1.5 

64 

547 

11° 

00 

21 

.0194 

.44 

65 

584 

18° 

01 

1 02 

: .0181 

.52 

77 

623 

0 

1 — I 

01 

53 

.0586 

.77 


value for autocorrelation analysis was allowed to be as high as one second. Yet 
at even this large lag value the computed autocorrelation showed a residual correla- 
tion of 10 to 20 percent. Note that the Fourier transform of a truncated function con- 
tains an additional (sin fA;f) term due to the discontinuity. The variable oj is the 
reciprocal of the lag value (transform variable) at the discontinuity. This explains the 
observed one to two Hz damped oscillations. Autocorrelation of actual data fell to 
zero typically within 0,8 seconds. Thus the residual correlation was apparently 
due to the theoretical analysis, probably the interferometer assumption. Figure 
54 is a plot of normalized autocorrelation computed for a single wind velocity, 

2 meter/sec, under the interferometer assumption. Note the abrupt switch to 
nearly constant value for Rj 9 ('T') at a lag value of 0.3 seconds. If contributions from 
all parts of the aperture were included instead of only two points at the extremes 
of the diameter, autocorrelation would be expected to go to zero due to an averag- 
ing effect. Part of the residual correlation may, however, be due to the 5/3 power 
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of radial separation in the phase structure function. This arises from the assump- 
tion of Obukhov -Kolmogorov turbulence spectra. For large separations, expressed 
in terms of lag value, and after conversion to image motion autocorrelation by 
division by k^(Vj^ r)^ the function decreases only as lag value to the negative one- 
third power. This results in a long tail to the theoretical autocorrelation function. 
Despite the oscillating behavior, the envelope of the fit spectral density followed 
the data well. For higher frequencies the oscillations tended to damp out. 

Turbulent profile parameters determined by the fit to image motion data were 
A, and a. Early in this process, a was fixed at 5/3. This gave reasonable 
results and also corresponded to the best fit to the Hufnagel curve. Ground 
strength of turbulence, A, scaled the computed spectral density to give the correct 
value of 0-2 . The exponential controlling parameter, H , determined power spectra 
shape. Typical values were in the range 100 meters to 1 kilometer. It was 
found that accuracy of fit to the data was not very sensitive to the value of Hq. 

The test statistic, sum of absolute differences, varied less than 10 percent over 
this range. This is in sharp contrast to irradiance differences which obtained a 
sharp minimum at the best fit values of B, , and/5 . Thus image motion spectral 
density, at least in this present form, is unable to give an accurate measure of 
the shape of lower altitude turbulence or else the form of the model does not lend 
itself to a sensitive test. The dominance of the low-frequency peak may explain 
this behavior. Uncertainties concerned with the height of this peak and its effect 
on model parameters are discussed in Appendix B. The present analysis does 
suggest that image motion effects result primarily from a region of high turbulence 
within a few hundred meters of the ground. In general, the best value of was 
larger for conditions of low wind near the surface, indicating a need to include 
higher winds available at the higher altitudes. Apparently they were needed to 
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construct the higher frequencies of image motion. Even so, turbulence beyond 
a few kilometers did not contribute sensibly to image motion. Although wind 
speeds at the tropopause showed a large variation from night to night, speeds for 
the first few kilometers were fairly uniform for most nights. This could explain 
why image motion spectral densities did not have significantly different shapes 
from night to night. The one exception is night 36 data. For this night, wind speeds 
as high as 30 meters per second were encountered in the first three kilometers 
with an unusual concentration of high winds in general at low altitudes. The 
corresponding image motion spectrum shows at least a factor of two increase in 
bandwidth and a lessening of the low-frequency peak. 

It is apparent that a profile could be chosen to generate data for good fits to 
all irradiance spectral densities. The computer generated fit, while determined 
primarily by the low-frequency data, followed the high-frequencies very well and 
tended to represent the average of the data points all along the curve. Little attempt 
was made by the computer fit to follow the data fine structure since there were 
not enough degrees of freedom in the model. Comparison of consecutive 20 
second records showed that the fine structure is for the most part a real effect and 
not statistical error. Note records one and two for data run 584 night 65 (Figures 
39 and 40). It is possible that equipment vibration or electronic noise could account 
for this fine structure, although the appearance of fine structure at different fre- 
quencies on different nights and different times argues against this. 

Perhaps the most striking feature is the maximum of irradiance spectral 
density that occurs between 10 and 20 Hz. In some cases, night 297 and night 65, 
the computer fit exhibits the same behavior. By printing out the tabulated 3 {• } 
function versus altitude for each frequency of interest, it was possible to determine 
which altitude ranges had the most weight in determining magnitude of spectral 
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density at a particular frequency. Low-frequency behavior, including that interest- 
ing maximum, was found to be controlled for all wind profiles primarily by con- 
tributions from the very high altitudes (above 15 kilometers). Contributions to 
the middle and higher frequencies were increasingly weighted to those altitudes 
where wind velocity maximized. Since low-frequency components were relatively 
stronger in terms of CIV^/Hz the turbulence profile was forced to maximize above 
the region of maximum wind. The maximum of C^(h) appeared at higher altitudes 
for nights with the largest maximum wind velocities (night 65 and 77) and the 
Gaussian "bump" was more sharply peaked. Low wind conditions represented by 
night 64 produced a broadly peaked C^(h) profile. This was to be expected since 
the lower wind was also more uniform with altitude. 
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CHAPTER VI 


CONCLUDING REMARKS 

The spectral density technique for estimating vertical profile of turbulence 
strength, Cj^(h), appears to be a significant improvement over previous methods 
based only on the magnitude of irradiance and image motion variance. The 
variance technique uses theoretical dependence on moments of the distribution. 

The spectral density technique also follows this approach with the moments being 
supplied by wind velocity as a function of altitude. A great deal more information 
is available, however, since wind velocity structure is not a simple power of alti- 
tude but provides a whole series of moments. Wind velocity values are accurately 
obtained from Weather Bureau rawinsonde flights adjacent to the site of stellar 
observations. These should be supplemented by local surface and near-surface 
measurements. By using near -zenith stellar observations the azimuth effects on 
wind velocity are eliminated. Thus wind velocity does not have to be a limiting 
factor in the determination of C^(h). The real limiting factors appear to be 
questions of statistical accuracy and lack of stationarity . In image motion there is 
an additional problem in separating instrumental from atmospheric effects. 
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In this report the theoretical base for stellar observations was reported and, 
where necessary, extended to account for operation of the Stellar Image Monitor 
apparatus. Certain major assumptions were involved in the theoretical base. 

(1) Rytov ’’approximation'' to the scalar wave equation. 

(2) Local stationarity and isotropy of optical statistics. 

(3) Obukhov - Kolmogorov turbulence spectra. 

(4) Taylor’s hypothesis regarding the relation of spatial and temporal fre- 
quencies through wind velocity. 

Assumptions (2) and (3) are the most contested. Nonstationarity effects were 
noted in the data and resulted in larger error bars as well as a divergence of the 
data from the theoretical description. Departures from an Obukov -Kolmogorov 
turbulence spectrum could also be expected for the stable near-surface conditions 
at night. 

Irradiance statistics seem to be well covered by the theory, although more 
attention is needed to the effects of aperture size and its interaction with zenith 
angle. Image motion statistics need more work, particularly in the description 
of image time autocorrelation function for finite apertures. It will be necessary 
to include phase effects from all regions of the given aperture. This will enable 
better computer fits to image motion spectral densities with a corresponding 
improvement in determination of profile. 

Data collection and reduction procedures should also be strengthened. In 
general the Stellar Image Monitor apparatus was found to be suitable for the types 
of data required in this thesis research. Irradiance statistics were measured with 
good accuracy and were available for the interesting range of aperture sizes near 
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the correlation distance. The use of larger diameter optics would provide new in- 
formation on image motion. In particular, a test of the weak aperture 

dependence from the Obukhov-Kolmogorov theory could be made. It would also 
be desirable to record and analyze longer records of image motion data. This 
would provide more information on low-frequency behavior, narrow the statistical 
error bars, and tend to average out effects of short-term nonstatlonarity in the 
atmosphere. Alternately, if data were taken with higher surface wind velocities, 
spectral densities would broaden and apparent resolution would improve. 

Computed C^(h) profiles as reported in the last chapter, although not incorporat- 
ing all the accuracy and improvements suggested, do provide valuable information 
and new insight into the shape of the actual profile. The following conclusions 
can be drawn from this work: 

(1) Nighttime contributions to (h) arise primarily from two turbulence re- 

N 

gions, one near the earth^s surface and another in the vicinity of the 
tropopause. 

(2) Image motion variance and spectral density can be explained adequately 
by contributions in the first few kilometers above the surface. If an 
exponential distribution for turbulence is used, the 1/e point is typically 
a few hundred meters. 

(3) The shape of irradiance spectral density and the structure of upper 
altitude winds require a definite contribution to C^(h) at altitudes near the 
tropopause. The assumption of a single maximum, Gaussian "bump", 
fits the data well for different wind profiles and total turbulent strengths. 
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(4) Parameters of the Gaussian fit tended to give a low, wide maximum 

for low wind conditions and a sharper, higher peak for high winds. This 
brings out the association of wind velocity and strength of turbulence 
implicit in the spectral density theory. Ts^pical maximum values of C^(h) 
near the tropopause were near 1 x lO m These were about 

three to four orders of magnitude less than the ground strength values, 
but supplied essentially all the contribution to irradiance statistics. 

It is apparent that statements in further detail cannot be made without in- 
creasing the model complexity and/or increasing data accuracy. To get most 
value out of refinements like this it would also be necessary to determine the 
appropriate number of degrees of freedom of the model, given data accuracy 
and given the method of model testing. 
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APPENDIX A 


IMAGE INTENSITY PROFILE AND OPTICAL TRANSFER FUNCTION 


The diffraction image formed in an optical system is characterized by a 
normalized intensity distribution called the optical spread function, 

S S(x, y; x' , y ^ = S(x, x^ \), 

light flux in + d\), which originates in dx' and passes through the optics to 
fall on dx, divided by dx dx' d\ , where 

X = two-dimensional vector in object plane 

and 

X = two-dimensional vector in image plane. 

The optical spread function is assumed to be time independent and normalized 
such that 



X 


; dx 


1 . 


Use of a transfer function to describe the optical system requires the exist- 
ence of an isoplanatism patch, an area in the object field over which wave aberra- 
tions are effectively constant. This enables the optical spread function to be 
written 

S - S(x“x ;X.)- (A2) 


The optical transfer function, T( f ) , is defined as the Fourier transform 
of S. The transform variable is f , two-dimensional spatial frequency, so that 
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T(f) 


I S(x - ; X.) exp [- 2t 7 i f * (x - x' )] dx 

m 


(A3) 


and 


S(x“x^;X) - j T(f) exp \^2u i f • (x-x')] df. 

- 00 


(A4) 


The object intensity distribution, 0(x' ) , can be written in terms of a Fourier 
integral of spatial frequencies^^ as 


I 


O(x') = W(f) exp (2 t 7 i f • x' ) df 


(A5) 


and 


W(f) 


^ r 


^ ) exp (- 2 t 7 i f • x' ) dx'. 


(A6) 


where 

W(f) = Fourier transform of O(x'). 

Image intensity distribution, I(x) , is by definition of the optical spread function, 
just the con solution of S(x - x'; \) with O(x'), that is, 


1 ( 5 ) 


^00 

I S(x~x'; /V)0(x')dx'. 


(A7) 


Equation (A4) for S and equation (A7) for I allow image intensity to be expressed 
in terms of the transfer function by 


1 ( 5 ) 



\ ) 0( X ' ) d5 ' 


(A8) 


00 

j j T( f ) exp [ 277 ! f- (x-x')]dfO(x')dx^ 
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00 

- ^ ^ T( f ) exp [ 277 - if • (x “ x' )] 0(x' ) dx' df 




T( f ) exp (2 tt i 




0( X ' ) exp (~ 277 i f • x')dx' df 




T( f ) exp (2 tt i f • x) W( f ) df . 


(A9) 


If the object is a point source, O(x') - S(x') and 


Thus, 


W(f) 


I(x) 


^00 

S(x 

- 00 


' ) exp ( - 2?/ i f • x' ) dx ' - 1 . 


^00 

- c 


T( f ) exp (27fi f • x)df 


(AlO) 


The variable, T(f), is, in general, complex and can be written in terms of a 
modulus and complex argument as 

T(f) = M(f) exp [in i c/)( f ) ] . 

For real symmetric images the transfer function is real, T( f ) = M( f ) , and 


r 


I(x) - M( f ) exp ( 2 t 7 i f • x) df 


(All) 


A transformation to polar coordinates for both x and f allows I to be expressed 
in a form more suitable for computation. That is, given 


00 

I(x, y) = JJ*M( f , g) exp [277 i ( fx + gy)] df dg , 
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and letting 
and 


then, 


i - CO cos <p , 
g = cosin(p. 


X - r cos 0 5 
y = r sin (9, 


^00 -277 

I( r, <9) = I I }A(^co, <p) exp [2tt i r co cos {(j) - 9) j co d(pdco . 

•'o •'o 


Symmetry arguments require that I be 9 independent and M be 0 independent. 
Therefore , 


poo p2?7 

I(r) = M( co) CO I exp ( 2t7 1 r oj cos 0) d<p dco , 

•^0 


or 


^ n 


I (r) = 2 77 I M(o;) o) Jq ( 2t 7 r o;) dco 
0 


(A12) 


This is the final form relating image intensity and the transfer function. 

The theoretical transfer function for a diffraction-limited optical system 
with circular aperture has been computed and results in 

_ 2 r - 1 \co [ ( kco\^ 

U(cv) - ^ icos 1 ^ 

where ^ 

X. = wavelength, 

D ■ aperture diameter , 

and 

CO = spatial frequency in cycles per arc second. 

Substitution of this expression into equation (A12) yields the familiar Airy pattern 
for a diffraction-limited lens: 



I(r) 



'2 


Jj (7rDr/\F)-|2 

ttPtTXF 
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where 


F = system focal length • 

In angular measure this becomes 


I(S) 



ttDS 

>v 


with S in radians. This expression is plotted in Figure 25. 


(A14) 
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APPENDIX B 


IMAGE MOTION LOW-FREQUENCY BEHAVIOR 

Image motion variance as measured in this research often exhibits larger 
variability among adjacent data runs than can be explained by electronic signal- 
to-noise ratios or statistical sample size error bars. In examination of image 
motion spectral density it is seen that the greatest variability as well as the 
greatest fraction of total variance (typically 30 to 50 percent) resides in the low 
frequencies below one Hz. In the consecutive records shown in Figures 39 and 40 
this low frequency peak is seen to change by a factor of four. Previously (Chapter 
IV) this has been attributed to nonstationarity in the atmosphere. This may very 
well be the case, for gusts of wind mixed with periods of calm could conceivably 
alter the total variance and low-frequency spectrum quite rapidly. In addition, a 
change in refraction contributed by the total atmospheric path would show up as a 
low-frequency effect. 

As resolution bandwidth of the spectral density analysis was increased from 
the nominal value of one Hz the large contribution at low frequencies was found 
to be concentrated in a smaller region, typically within 0.3 Hz of the origin, 
and to increase in height. Contributions in this region may be due entirely to an 
instrument effect, tracking of the star by the telescope drive. The standard 
computer analysis described in Chapter IV removed any linear component in the 
image motion data by least squares fitting a straight line to the data. This was 
done to account for drift in the tracking. It is also possible that improper tracking 
could contribute curvature to image motion data. This very point is discussed 
by Kolchinskii . He concludes that, while small, the curvature component could be 
significant for stars away from the meridian and in the case of image wander 
caused by changes in total refraction. 
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To test this assertion, a quadratic term was also removed in a least squares 
fashion from a sample of data showing a strong low-frequency peak, data run 
558 of night 64. This procedure reduced the peak a factor of two and resulted 
in a more believable shape for spectral density, i.e. , a gentle rolling off as zero 
frequency was approached. Removal of higher order terms did not further signif- 
icantly reduce the peak, indicating that only the linear and quadratic terms were 
appreciable in the data. It is still difficult to state how much of this effect is 
truly an instrument effect. The data could very well contain atmospheric contri- 
butions at these low frequencies that would be attenuated by removal of even the 
linear term. In fact, the theoretical prediction of the Obukhov -Kolmogorov spectrum 
is that for low winds (less than one meter per second) significant correlation 
exists to lag values beyond one second of time. This indicates significant contri- 
bution at frequencies below one Hz. Evidence of this is supplied in the relative 
ease with which a computed spectral density was constructed for a reasonable 
C^(h) profile to fit the low-frequency peak. 

For present purposes suffice it to say that for resolution bandwidth of one 
Hz the spectral density values at zero and 0.5 Hz should be considered accurate 
only to within a factor of 5 to 10, and are probably a little high in general. The 
net effect of higher than actual peak on model parameters for Cj^(h) is to require 
a smaller H^ and larger A. This is equivalent to concentrating turbulence at 
lower altitudes where wind speed is lower and low-frequency image motion can 
be enhanced. Add to this the large uncertainty and variability of low altitude 
winds and it is difficult to make more than order of magnitude estimates of the low- 
altitude turbulence distribution. 



APPENDIX C 


IMAGE SIZE CALCULATION 

The image size determination of equation (75) is greatly simplified mathe- 
matically if the approximation of exp (- cf ^ is made for the true exp (- ) 

MTF factor. To test this approximation a numerical integration of equation (75) 
was carried out for both cases. This integration was performed for one value 
of the parameter C which generated an MTF close to that of stellar data. With 
f expressed in cycles per arc second the value used for C was two. Numerical 
integration was performed on a digital computer through a standard integration 
subroutine. The subroutine employed the trapezoidal rule together with Rom- 
berg’s extrapolation method. Although the formal integration extended from 
zero to infinity the MTF exponential factor in the integrand enabled the numeri- 
cal procedure to be truncated at a reasonable value for spatial frequency. Thus 
the actual integration was carried out over the range zero to 50 cycles per arc 
second. 

The integrals were performed successively for different values of the vari- 
able S, arc seconds, and image profiles were constructed for both cases. The 
square law approximation resulted in a slightly wider image with less intensity 
in the tail. Despite this behavior the two curves followed each other quite closely 
and the total error in rms image width was only 7 percent, with the square 
law case giving the higher value for width. Thus it appears reasonable to employ 
the square law approximation. 
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